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Abstract 

When  dealing  with  acoustic  signals,  the  environment  in  which  the  signal  propagates 
affects  the  signal  in  measurable  ways.  These  effects  lead  to  echoes,  changes  in  amplitude, 
and  sometimes  a  confluence  of  other  signals  which  may  render  the  signal  useless  in  terms  of 
information  retrieval.  When  considering  acoustics  in  terms  of  pressure  changes  due  to  the 
driving  signal,  we  are  able  to  model  and  measure  the  effects  of  the  environment  on  a  given 
signal.  Once  this  model  is  known,  we  are  able  to  completely  invert  the  environmental 
transfer  enacted  on  the  signal  and  filter  out  other  signals,  as  long  as  some  assumptions  are 
held  in  the  implementation  of  this  procedure. 
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ENVIRONMENTAL  ACOUSTIC  TRANSFER  FUNCTIONS  AND  THE 


FILTERING  OF  ACOUSTIC  SIGNALS 

I.  Introduction 

1.1  Background 

Signal  processing  is  the  method  of  taking  a  given  signal  and  extracting  useful  in¬ 
formation,  usually  by  a  means  of  a  transformation  of  some  kind.  Acoustic  signals  are 
functions  of  time  in  which  the  output  is  a  pressure  or  a  velocity  potential  response.  An 
acoustic  signal  is  affected  by  the  environment  in  which  it  propagates,  so  one  can  attempt  to 
remove  the  environmental  effects  to  extract  the  useful  information,  in  this  case  the  original 
signal.  The  classic  example  of  this  problem  and  how  the  human  brain  deals  with  it  is 
called  the  Cocktail  Party  Problem  (3).  Imagine  yourself  in  a  room  full  of  people,  where 
different  conversations  are  happening  all  around  you,  and  you  attempt  to  either  participate 
in  a  conversation  around  you  or  listen  to  a  conversation  near  you.  Somehow,  your  brain 
processes  this  multitude  of  signals  by  localization  techniques  and  filters  the  information 
in  such  a  way  as  you  are  able  to  ignore  these  extraneous  conversations  and  concentrate 
on  a  particular  conversation.  This  thesis  will  derive,  in  a  mathematical  framework,  the 
process  of  filtering  extraneous  signals  in  a  way  that  yields  the  original  signal,  and  will  then 
apply  this  process  to  the  Cocktail  Party  Problem  in  an  attempt  to  describe  how  useful  this 
ability  can  be. 

This  ability  has  many  applications  in  both  the  Department  of  Defense  and  commer¬ 
cial  industries.  Considering  an  acoustic  signal  to  be  a  form  of  information,  the  basic 
purpose  of  filtering  the  received  signal  is  to  retrieve  the  original  information,  such  as  in 
intelligence  gathering.  One  use  for  the  Department  of  Defense  is  the  use  in  "bugging"  a 
room  with  listening  devices  that  record  the  information  being  disseminated  in,  for  example, 
a  terrorist  meeting.  If  a  number  of  people  conversing  in  a  meeting,  the  received  signal 
will  be  the  sum  of  these  signals  at  their  respective  locations.  Taking  into  consideration 
their  different  locations,  the  recording  can  begin  to  sound  quite  garbled  and  unintelligible. 
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Since  this  process  is  being  treated  as  information  retrieval,  this  garbled  recording  would 
simply  be  useless  information.  By  finding  a  method  for  separating  the  signals  and  remov¬ 
ing  environmental  effects,  the  original  information  would  then  be  recovered,  potentially 
containing  useful  information. 

The  private  sector  also  has  many  uses  for  this  method,  and  while  none  of  them 
may  not  have  as  direct  an  impact  on  our  daily  lives  as  the  defense  application,  they  can 
be  helpful  at  a  relatively  low  cost  from  a  computing  standpoint.  Consider  a  boardroom 
meeting  transpiring  in  a  room.  While  at  most  times,  there  will  only  be  one  person  talking, 
there  may  be  at  other  times  many  people  talking  at  once,  arguing  with  each  other.  This  is 
when  the  application  of  this  signal  processing  will  become  most  useful.  The  entire  meeting 
could  be  recorded  with  every  word  clearly  spoken  and  usable  for  future  reference.  Also, 
consider  the  explosive  proliferation  of  cellular  phones  in  the  past  decade.  While  in  most 
cases,  cellular  phones  are  easy  to  use,  there  are  some  noisy  environments  in  which  they  are 
practically  useless.  The  user  of  the  phone  can  eliminate  this  problem  by  using  an  earpiece, 
sending  the  signal  directly  to  his  or  her  ear,  but  the  person  on  the  other  end  of  the  line  is 
still  stuck  with  the  multitude  of  signals  entering  the  mouthpiece  of  the  phone.  One  way 
to  eliminate  these  sounds  would  be  to  filter  the  other  voices  out  of  the  incoming  signal,  by 
somehow  correlating  the  user’s  voice  with  an  amplitude  function  or  some  other  way  and 
to  then  filter  the  sounds  by  software  included  in  the  phone.  One  product  currently  on 
the  market  that  performs  this  task  is  a  mobile  phone  headset  called  Jawbone.  Jawbone 
uses  technology  developed  by  DARPA  to  correlate  the  motion  of  the  user’s  jawbone  with 
the  signal  being  received.  When  the  user’s  jawbone  is  not  moving,  the  headset  effectively 
turns  off  and  uses  the  ambient  sound  being  received  during  this  time  to  create  a  noise 
filter  for  use  when  the  user  is  talking.  This  is  a  perfect  example  of  this  technology  being 
marketed  in  the  private  sector  (1). 

1.2  Definitions 

This  section  presents  and  defines  certain  terms  which  are  used  throughout  this  thesis. 
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Environment  refers  to  the  finite  region  of  influence  acting  on  a  signal.  In  the  math¬ 
ematical  depiction  of  the  problem,  this  will  be  replaced  by  the  domain  with  boundary 
conditions. 

Attenuation  is  the  effective  reduction  in  amplitude  of  a  signal  caused  by  environ¬ 
mental  factors.  For  the  purposes  of  this  thesis,  attenuation  will  be  described  as  either 
atmospheric  or  boundary. 

Isothermic  describes  the  property  of  a  medium  having  a  uniform  temperature,  i.e. 
the  temperature  gradient  with  respect  to  the  three  spatial  variables  and  time  is  equal  to 
zero. 

Isobaric  describes  the  property  of  a  medium  having  a  uniform  pressure  density,  i.e. 
the  dispersion  of  the  medium’s  constituent  gases  is  constant  throughout  the  environment. 

The  coefficient  of  reflection  is  the  ratio  of  the  energy  of  a  wave  after  encountering  a 
barrier  to  its  original  energy.  While  the  coefficient  of  reflection  can  take  complex  values, 
the  examples  in  this  thesis  deal  only  with  real-valued  coefficients.  This  implies  that 
there  is  no  phase  shift,  or  modulation,  on  the  signal  when  encountering  a  boundary,  only 
attenuation. 

Ultrasonic  describes  a  frequency  above  the  human  threshold  of  hearing. 

A  Source  is  the  origin  of  a  signal,  which  will  be  a  person  speaking,  unless  stated 
otherwise. 

A  Receiver  is  the  instrument  with  which  a  recording  of  an  environmentally  altered 
signal  will  be  made.  For  the  purposes  of  this  thesis,  a  receiver  will  refer  to  an  omnidirec¬ 
tional  microphone  unless  stated  otherwise. 

Transfer  Function  refers  to  the  rule  by  which  an  environment  affects  a  signal.  It  is 
our  goal  to  find  the  transfer  function  of  a  given  environment. 

1.3  The  Problem 

An  acoustic  signal  can  be  distorted  in  a  number  of  ways  by  an  environment,  thereby 
making  information  contained  within  the  signal  less  accessible.  The  effect  of  the  environ- 
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ment  must  be  determined  in  order  to  retrieve  the  original  information.  This  thesis  pursues 
those  goals  by  finding  solutions  to  the  following  three  problems: 

1.  Problem  1  is  defined  as  that  of  finding  and  inverting  the  room  transfer  function. 

2.  Problem  2  is  defined  as  the  problem  of  filtering  different  signals  being  transmitted 

at  once. 

3.  Problem  3  is  defined  as  retrieving  the  signal  of  a  mobile  source. 

1.4  Research  Objectives 

This  thesis  develops  a  rule  for  environmental  effects  on  a  signal,  as  well  as  a  rule 
for  removing  these  effects.  This  removing,  or  filtering,  enables  us,  to  a  certain  degree,  to 
separate  out  multiple  time-dependent  signals  into  separate,  whole  signals  containing  the 
original  information  contained  in  each  signal. 

1.5  Scope 

The  results  of  this  thesis  enable  us  to  determine  under  which  conditions  a  filtering 
operation  can  successfully  be  performed  on  a  set  of  received  signals.  While  Problem  2  is 
primarily  concerned  with  the  two-receiver /two-source  case,  this  can  be  applied  to  numerous 
other  cases,  such  as  multiple  (more  than  two)  signals  and  signals  with  a  variable  location 
(from  here  on,  referred  to  as  mobile  sources  for  the  purpose  of  this  thesis). 

The  results  of  Problem  2  also  lead  to  a  development  of  a  rule  for  judging  how  "good" 
a  filter  is,  based  on  the  condition  number  of  the  matrix  of  transfer  functions. 

While  the  scope  of  the  thesis  seems  somewhat  narrow  at  the  beginning,  the  findings 
can  be  broadened  to  apply  to  a  general  case  (the  most  general  of  which  is  the  Cocktail  Party 
Problem  with  no  limitations  on  source  locations).  The  benefit  of  directional  microphones 
in  this  setting  is  discussed,  but  not  implemented  in  the  mathematics. 
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1.6  Approach 

The  approach  to  the  very  basic  problem  of  finding  the  room  transfer  function  is 
that  of  solving  the  three-dimensional  wave  equation  modified  to  incorporate  atmospheric 
attenuation.  The  boundary  conditions  and  initial  conditions  are  developed  in  the  latter 
part  of  Chapter  II. 

1.7  Thesis  Outline 

This  thesis  takes  a  bottom-up  approach.  Chapter  II  deals  with  a  literary  review  of 
most  of  the  concepts  and  topics  discussed  in  this  introductory  chapter.  Chapter  III  begins 
with  a  statement  of  the  foundational  assumptions  made  in  solving  this  problem,  as  well 
as  detailing  the  solution  to  the  initial  problem.  Finally,  the  results  of  the  core  problem 
of  finding  the  room  transfer  function  are  extended  in  the  same  mathematical  framework 
to  the  filtering  problem  (Problem  2)  and  the  tracking  problem  (Problem  3).  Chapter  IV 
discusses  the  implementation  of  this  process  on  a  specific  situation  where  numerical  values 
are  assigned  to  the  parameters  used  in  the  derivation  of  the  solution  to  the  initial  problem. 
Chapter  V  summarizes  the  findings  of  this  thesis  as  well  as  presenting  recommendations 
for  future  research  in  this  area. 
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II.  Literature  Review 


2.1  Introduction 

This  chapter  consists  of  a  literary  review  of  the  two  basic  aspects  of  the  three  prob¬ 
lems.  The  first  part  consists  of  a  review  of  the  fundamental  aspects  of  acoustics  and  how 
different  environments  affect  a  signal.  This  will  deal  primarily  with  part  of  Problem  1, 
that  of  finding  the  room  transfer  function.  The  second  part  will  consist  of  a  review  of 
filtering  techniques  and  processes,  which  generalizes  the  transformation  used  in  Problem 
1  to  a  more  complicated  system. 

2.2  Fundamentals  of  Acoustics 

2.2.1  A  Definition  of  Sound.  Before  beginning  the  goal  of  determining  a  rule 
by  which  an  environment  affects  an  acoustic  signal,  a  description  of  sound  must  be  given. 
"The  whole  study  of  sound  is  the  study  of  vibrations"  (8:36).  By  this  one  sentence,  written 
in  one  of  the  most  prolific  texts  on  acoustics,  a  great  deal  can  be  deduced.  The  key 
word  in  this  sentence  is  "vibrations".  When  considering  the  mathematical  description  of 
vibrations,  the  likely  tool  for  describing  them  will  be  the  trigonometric  functions  because  of 
their  oscillating  qualities.  Morse  begins  his  description  of  what  sound  really  is  by  speaking 
of  an  oscillating  system,  using  the  example  of  a  spring-mass  oscillator.  By  introducing  this 
example,  he  shows  that  in  a  real-world  oscillating  system,  there  is  dampening  caused  by 
friction  (8:37).  In  an  acoustic  environment,  this  dampening  is  the  atmospheric  attenuation. 

2.2.2  Measuring  Sounds.  Now  that  sound  has  been  defined  the  system  in  which 
it  propagates  has  been  described,  a  method  to  quantify  these  effects  must  be  discussed. 
Thus  far,  the  environment  has  simply  been  described  as  air  and  been  called  a  medium  of 
propagation.  Since  air  is  composed  of  gases,  it  follows  the  laws  of  fluid  dynamics,  and 
so  these  laws  must  be  used  to  describe  the  effect  of  a  signal  on  the  air.  A  fluid  has  a 
number  of  characteristics  that  determine  how  it  behaves.  In  this  case,  the  fluid  density, 
the  pressure,  and  the  temperature  are  the  key  factors  that  determine  how  a  wave  will 
propagate  through  this  fluid.  When  discussing  the  propagation  of  a  wave,  there  must  be 
a  way  to  quantify  the  magnitude  of  the  wave’s  effects  at  a  spatial  point  in  the  environment 
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at  a  certain  time.  Since  the  direct  measurement  of  temperature  and  density  is  rather 
difficult,  it  is  more  efficient  to  measure  the  amplitude  by  its  pressure.  Measurements 
of  pressure  are  also  more  accurate  than  those  of  temperature  and  density  (8:229).  The 
familiar  measurement  of  the  intensity  of  sound,  decibels,  is  a  logarithmic  measurement  of 
pressure.  It  is  in  this  framework  of  pressure  changes  caused  by  sound  waves  that  a  signal 
and  the  environment’s  effects  are  measured. 

2.2.3  Deriving  the  Wave  Equation.  Speaking  in  the  framework  of  pressure, 
Morse  derives  the  wave  equation  from  a  succession  of  equations  used  to  describe  the  law 
of  conservation  of  mass  and  the  law  of  conservation  of  momentum  in  a  fluid.  For  a  fluid, 
the  law  of  conservation  of  mass  in  three  spatial  dimensions  is 

^  +  V-pu  =  0  (2.1) 

and  the  law  of  conservation  of  momentum  is 

+  Vp  =  0,  (2.2) 

where  p  =  p(x,  t )  is  the  particle  density  of  the  fluid  at  location  x  at  time  t,  V  (x,  t )  is 
the  particle  velocity  in  the  fluid,  and  p  (. x ,  t)  is  the  fluid  pressure.  Since  the  environment 
consists  of  a  homogeneous  fluid,  then 


dV 

dt 


+  {p-  V)v 


p(x,t)=p  =  c  (2.3) 

for  some  constant  c  >  0.  In  a  fluid,  the  pressure  is  a  function  of  the  particle  density,  so 
define  /  such  that 

P  =  f(p)  •  (2-4) 
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This  definition  is  called  the  equation  of  state.  These  parameters  can  be  linearized  about 
initial  conditions  such  that 


V  (x,  t ) 

=  Po+P  (x,  t ) 

(2.5) 

P 

=  Po  +  P 

(2.6) 

V  ( x ,  t ) 

=  Vq  +  P  (x,  t ) 

(2.7) 

where  vq  =  0  since  the  initial  conditions  are  quiescent.  Using  these  linearized  functions, 
the  equations  of  mass,  momentum,  and  state  can  be  restated  as 


=  o 

_  „  di> 

Vp+p°at  =  0 

p  =  C2p  . 


(2.8) 

(2.9) 

(2.10) 


Taking  the  derivative  with  respect  to  time  of  the  equation  of  mass  above,  the  following 
equations  are  found 


d2p  d 

W  +  m  (p°v ' v) 

d2p  d 

a?  “  at («>v ' v) 


0 

Ap  , 


(2.11) 

(2.12) 


and  the  equation  relating  the  particle  density  and  pressure  emerges  as 


d2  p 

w  =  Ap 


(2.13) 


but  the  equation  of  state  can  be  resubstituted  back  into  this  equation,  leaving  the  familiar 
wave  equation 

A p  -  4 ptt  =  0  .  (2.14) 

cz 

2.2.4  Attenuation- Atmospheric  and  Boundary.  In  Chapter  I,  attenuation  is 
defined  as  the  reduction  in  amplitude  of  an  acoustic  signal.  Going  back  to  Morse’s  example 
of  an  oscillating  spring,  he  describes  the  friction  in  this  system  as  a  resistance  to  a  change 
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in  inertia  of  the  medium.  In  the  acoustic  setting,  this  friction  is  the  resistance  of  air 
molecules  to  a  displacement,  thereby  reducing  the  energy  of  the  signal  by  a  reduction  in 
amplitude,  which  is  the  definition  of  attenuation.  The  amount  of  attenuation  is  given  by 
a  vector  T,  which  defines  the  amount  of  attenuation  in  the  directions  of  the  domain.  In 
this  case,  the  vector  has  constant  components,  all  of  which  are  equal  such  that 


r  =  2T 


(2.15) 


This  implies  that  there  is  no  directional  dependence  in  atmospheric  attenuation.  This 
assumption  is  justified  in  that  the  propagational  environment  is  isothermal  and  isobaric, 
with  no  fluid  flow.  Since  this  attenuation  takes  places  as  the  signal  propagates  through 
space,  its  effect  on  the  pressure  in  the  sound  wave  can  be  described  using  the  gradient  of 
the  pressure.  The  pressure  caused  by  the  wave  is  decreased  as  the  wave  passes  through 
space,  and  so  the  atmospheric  attenuation  can  be  given  by  the  expression 


f  •  Vp  =  2T  (px  +py+  pz ) . 


(2.16) 


In  the  differential  equation  developed  in  later  sections,  this  term  is  added  to  the  wave  equa¬ 
tion,  giving  a  modified  form  for  the  wave  equation  which  takes  into  account  atmospheric 
attenuation. 

The  atmospheric  attenuation  is  not  the  only  source  of  attenuation,  however-in  this 
setting  there  also  exists  boundary  attenuation.  In  the  case  of  the  oscillating  spring  mass 
system,  consider  a  second  spring  set  at  the  lowest  point  at  which  the  oscillating  spring 
will  reach  in  an  oscillation.  When  the  mass  interacts  with  the  boundary,  it  will  encounter 
resistance  to  any  further  displacement  by  this  second  spring.  This  resistance  will  effectively 
reduce  the  energy  of  the  oscillating  spring,  thereby  reducing  the  amplitude  of  the  next 
oscillation.  Now,  translating  this  idea  to  the  acoustic  setting,  the  signal  will  reach  a 
boundary  and  will  find  a  resistance  to  a  change  in  inertia.  This  resistance  will,  analogous 
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to  the  oscillating  system,  cause  a  reduction  in  amplitude,  thereby  attenuating  the  signal 
at  the  boundary. 

The  comparison  to  a  spring  on  the  boundary  is  fairly  accurate.  The  spring  constant 
translates  to  an  impedance  on  the  boundary,  which  characterizes  how  much  sound  the 
boundary  absorbs  and  how  much  is  reflected.  If  we  assume  that  the  boundary  reflects 
a  signal  with  constant  attenuation,  that  is,  no  frequency  dependence  and  no  phase  shift, 
then  the  coefficient  of  reffection  of  a  boundary  is  a  real,  positive  constant.  A  reflective 
coefficient  of  (3  leads  to  the  Robin  boundary  condition 

£=»  <*■”> 

where  denotes  the  outward  unit  normal  derivative  of  pressure  with  respect  to  the 
boundary.  The  left-hand  side  is  the  reflected  signal  and  the  right-hand  side  is  the  original 
signal. 

A  coefficient  of  reffection  approaching  infinity  implies  that  the  boundary  absorbs 
all  of  the  incident  wave  energy  and  none  is  reflected.  This  complete  attenuation  of  the 
incident  signal  is  described  by 

p  =  0.  (2.18) 

The  most  realistic  case  is  that  of  a  boundary  that  absorbs  some,  but  not  all  of  the 
signal,  and  so  the  reflected  wave  undergoes  a  reduction  in  amplitude.  This  boundary 
condition  is  given  by 

|  =  ft>  (2.19) 

where  /3  is  the  reflective  coefficient  of  the  boundary.  The  reffection  coefficient  can  be  a 
function  of  frequency,  attenuating  different  frequencies  by  different  amounts,  but  for  the 
purposes  of  this  thesis,  j3  is  a  positive,  real  constant. 

2.2.5  Forcing  Functions.  Consider  an  atmospheric  environment  with  isothermal 
and  isobaric  properties.  Since  the  temperature  and  pressure  gradients  are  both  equal  to 
zero,  then  the  system  will  be  at  rest  until  another  force  acts  on  the  environment  if  the  initial 
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conditions  are  quiescent.  In  this  case,  there  is  no  sound  being  propagated.  However,  if  the 
environment  is  acted  upon  by  another  force,  an  acoustic  signal,  then  the  environment  will 
be  disturbed  from  its  quiescent  state  and  a  signal  will  be  propagated  through  the  medium. 
This  force  can  be  considered  mathematically  as  a  driving  function  to  the  system. 

The  shape  of  the  source  partially  determines  the  nature  of  the  sound  wave  being 
propagated.  In  the  case  of  a  point  source,  the  sound  wave  propagated  is  spherical  in 
geometry  since  there  is  no  apparatus  to  act  as  a  focusing  antenna,  as  the  human  mouth 
does,  but  for  the  sake  of  simplicity,  this  thesis  will  deal  with  point  sources  that  will  radiate 
spherically  from  the  source  location.  This  point  source  can  be  described  as  a  time  depen¬ 
dent  signal  being  propagated  at  a  location  xq,  and  so  the  use  of  the  Dirac  delta  function 
is  appropriate  in  the  use  of  a  forcing  function.  A  time-dependent  signal  propagated  from 
a  location  xq  therefore  acts  as  a  forcing  function  in  the  form 

s(t)5(x  —  x o) .  (2.20) 

2.2.6  The  Modified  Wave  Equation.  Now  that  the  different  effects  of  the  envi¬ 
ronment  and  the  nature  of  the  source  have  been  discussed,  it  is  possible  to  formulate  a 
differential  equation  that  combines  these  effects.  It  has  been  mentioned  that  the  effect  of 
atmospheric  attenuation  must  be  included  in  the  differential  operator,  and  so  the  left-hand 
side  of  the  differential  equation  takes  on  the  form 

Ap  +  f-Vp-  4 ptt.  (2.21) 

The  forcing  function  makes  up  the  right-hand  side  of  the  equation,  and  so  the  differential 
equation  is 

Ap  +  f  •  Vp  -  4 Ptt  =  s(t)5(x  —  x o) ,  feii,  t  >  0  (2.22) 

cz 

The  effects  on  the  signal  at  the  boundary  must  satisfy  the  Robin  boundary  condition 

=  f3p,  x  €  <9D,  t  >  0.  (2.23) 
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and  must  have  initial  conditions 


p(x,t )  =  0,  x  £  fl,  t  =  0  (2.24) 

Pt  (x,t)  =  0,  x  £  fl,  t  =  0. 

These  quiescent  initial  conditions  are  justified  by  the  fact  that  the  microphones  can  be 

calibrated  such  that  the  pre-existing  pressure  in  the  room  is  ignored,  and  so  only  the 

changes  in  pressure  are  measured. 

2.3  Inversion  and  Filtering 

2.3.1  Localization  Cues  and  Qualitative  Elements  of  a  Signal.  When  a  signal  is 
acted  upon  by  the  environment,  there  are  many  changes  that  occur  in  the  original  signal, 
all  of  which  are  described  in  the  differential  operator  and  boundary  conditions.  These 
effects  reflect  the  nature  of  the  propagational  environment  and  are  instrumental  in  the 
location  of  the  source’s  origin.  Humans  have  auditory  receivers  in  the  form  of  ears.  The 
size  and  shape  of  these  ears  help  to  directionalize  the  signal  in  that  they  are  not  open 
all  around,  but  have  a  bowl  shape  to  them.  The  cartilage-composed  outer  part  of  the 
ear  is  called  pinna,  and  they  help  humans  to  determine  where  a  signal  is  coming  from-its 
azimuth  and  elevation.  Since  this  model  assumes  that  the  microphones  being  used  are 
omnidirectional,  the  azimuth  and  elevation  of  the  source  signal  is  not  taken  into  account, 
and  so  these  effects  on  the  signal  are  not  measured.  Only  the  source  and  receiver  locations 
with  respect  to  boundaries  and  the  distance  between  source  and  receiver  affect  the  signal. 

The  nature  of  the  environmental  boundaries  play  a  large  role  in  the  effect  of  the 
environment.  If  these  boundaries  are  highly  reflective,  the  signal  will  not  be  attenuated 
as  much  by  the  boundaries  and  so  atmospheric  attenuation  will  play  a  larger  role.  Since 
atmospheric  attenuation  generally  has  a  much  smaller  effect  than  boundary  attenuation, 
the  received  signal  typically  has  fewer  echoes,  but  if  the  walls  are  highly  reflective,  echoes 
will  be  more  prominent,  making  localization  more  difficult  (9). 
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2.3.2  Source  Localization.  One  of  the  main  assumptions  made  in  this  thesis 
is  that  the  receiver  and  source  locations  are  known.  In  the  boardroom  example,  this 
assumption  is  usually  valid,  but  in  the  intelligence  gathering  example,  the  source  locations 
will  be  harder  to  determine.  The  process  of  determining  source  positions  is  called  source 
localization. 

There  are  two  types  of  source  localization-active  localization  and  passive  localization. 
Active  localization  introduces  some  form  of  energy  into  the  system  to  gather  information. 
A  familiar  form  of  active  localization  is  found  in  remote  sensing-the  radar.  The  radar  is  the 
electromagnetic  counterpart  to  sonar,  which  uses  ultrasonic  pulses  that  feed  information 
about  the  environment  back  to  a  receiver.  Ultrasonic  imaging  is  widely  used  during 
pregnancies  to  observe  a  fetus  during  gestation.  It  can  also  be  used  in  the  gaseous 
environment  of  the  atmosphere,  but  since  the  atmosphere  is  much  more  tenuous  than  the 
human  body,  it  can  be  harder  to  focus  the  image  (2). 

An  inherent  problem  in  active  sensing  is  the  possibility  of  detection.  Since  an  active 
sensing  method  introduces  excess  energy  into  a  system,  this  excess  energy,  in  whatever 
form  it  may  be,  may  be  detected,  and  so  a  more  passive  method  of  determining  source 
localization  is  desired. 

If  there  are  two  microphones  in  a  room,  a  signal  will  take  different  lengths  of  time 
to  propagate  to  microphones  at  different  distances  from  the  source.  This  duration  is 
called  the  time  difference  of  arrival  (TDOA),  and  it  can  be  exploited  to  determine  source 
locations.  The  set  of  possible  locations  of  a  source  is  a  hyperbola  with  its  receiver  as 
the  focus.  Using  the  TDOA  method,  the  intersections  of  these  hyperbolae  can  be  found, 
which  will  give  the  location  of  the  source.  Since  this  method  uses  the  recording  as  the 
means  of  determining  source  locations,  there  is  no  further  preparatory  effort  needed.  If 
the  receiver  is  not  detected,  then  neither  is  the  process  of  source  localization  (5). 

2.3.3  Time-Reversal.  Once  the  source  location  has  been  determined,  the  informa¬ 
tion  needed  to  make  the  necessary  transformation  to  remove  the  effects  of  the  environment 
is  complete.  One  method  similar  to  the  attempts  being  made  in  this  thesis  is  called  Time 
Reversal.  This  is  a  method  for  removing  the  environmental  effects  of  the  propagational 
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medium  in  ultrasonic  images,  but  the  process  is  very  similar  to  that  of  filtering  an  acoustic 
signal.  When  ultrasonic  images  are  made,  the  environment  changes  the  reflected  signals 
making  the  images  blurred.  The  objective  is  to  remove  the  environmental  effects  on  the 
image.  When  this  process  is  put  into  the  framework  of  information  gathering,  time  re¬ 
versal  and  RTF  inversion  are  virtually  equivalent.  The  information  is  acted  upon  by  the 
environment,  giving  a  received  signal  to  which  a  transformation  is  applied  to  retrieve  the 
original  information.  The  time  reversed  focused  image  is  simply  the  filtered  signal-only 
the  nature  of  the  information  is  different  (6). 
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III.  Mathematical  Analysis 


3.1  Introduction 

In  this  chapter,  the  mathematics  of  acoustics  lead  to  solutions  for  the  three  problems 
posed  in  Chapter  I.  Initially,  a  general  solution  for  an  arbitrary  domain  will  be  explored 
and  then  the  spatial  operator  specified  to  a  rectangular  domain.  The  second  part  of  this 
chapter  deals  with  the  dimensional  analysis  of  the  wave  equation  in  an  attempt  to  ascertain 
exactly  what  physical  quantity  is  being  measured.  Section  four  consists  of  restating  of 
the  problem  by  Duhamel’s  Principle  so  that  the  differential  equation  can  be  solved  more 
readily  by  analytic  methods.  Sections  five,  six,  and  seven  deal  with  Problems  1,  2,  and 
3,  respectively.  In  each  of  these  sections,  it  is  shown  how  the  three-dimensional  modified 
wave  equation  is  the  key  to  solving  the  three  problems,  if  applied  properly. 

3.2  The  Eigenfunction  Approach 

Consider  an  arbitrary  domain  in  which  the  sound  wave  will  propagate.  Associated 
with  this  domain  are  the  boundary  conditions  that  define  the  behavior  of  the  wave  at  the 
boundary  of  this  domain.  Let  L  be  the  linear  differential  operator  that  will  define  the 
spatial  shape  of  the  wave  modes.  Then  L  is  the  spatial  operator  of  the  wave  equation  and 
so  the  homogeneous  wave  equation  in  this  domain  (scaled  with  respect  to  the  speed  of  the 
wave,  c)  can  be  written  as 

LU  =  Utt.  (3.1) 

where  U  =  U  (x,  t )  is  the  displacement  of  whatever  measurement  is  being  taken  -in  terms 
of  sound,  this  is  usually  pressure.  Assuming  that  the  operator  L  has  no  dependency  on 
t.  such  as  time-dependent  coefficients,  then  U  (x,  t )  can  be  assumed  to  be  separable  such 
that 

U  (. x ,  t )  =  u  (x)  T  ( t )  . 

Even  though  Chapter  II  defined  the  forcing  function  driving  the  wave  equation,  it  is  omitted 
here  because  of  the  application  of  Duhamel’s  Principle  in  section  four,  which  reformulates 
the  inhomogeneous  initial  boundary  value  problem  with  quiescent  initial  conditions  into  one 
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that  is  homogeneous  with  nonquiescent  initial  conditions.  By  finding  the  eigenfunctions 
un  of  the  spatial  linear  operator  L,  then  the  substitution 

LUn  =  -X nun  =  -A nun  (x)  T  (t)  (3.2) 

where  \n  is  the  eigenvalue  associated  with  the  eigenfunction  un,  eliminates  the  operator 
in  the  wave  equation,  yielding 


-A nunT  ( t )  =  unT"  (t)  .  (3.3) 

The  eigenvalues  are  all  positive,  and  so  therefore  there  is  a  nontrivial  solution  to  this  equa¬ 
tion  (7:389).  Dividing  both  sides  of  the  equation  by  un  eliminates  the  spatial  dependence 
and  yields  only  a  time-dependent  ordinary  differential  equation 

T"  +  k2nT  =  0  (3.4) 

k2n  =  A  n,  (3.5) 


which  has  the  solution 


Tn  ( t )  =  An  sin  (knt)  +  Bn  cos  (knt)  .  (3.6) 

The  solution  to  the  wave  equation  in  terms  of  the  eigenfunctions  of  the  spatial  operator 
L  is  therefore  the  superposition  of  the  individual  eigenfunctions  multiplied  by  the  time- 
dependent  solution: 

Un  (. X ,  t)  =  un  (x)  Tn  (t)  =  un  (x)  [. An  sin  (knt)  +  Bn  cos  (knt)]  (3.7) 

U  (x,  t.)  =  y,  Un  (x,  t)  .  (3.8) 

n 

From  this  point,  the  problem  becomes  that  of  finding  the  eigenfunctions  of  the  spatial 
operator  L.  The  operator  and  the  associated  boundary  conditions  determine  both  the 
shape  the  eigenfunctions  take  as  well  as  the  values  which  the  eigenvalues  An  can  take  on. 


3-2 


The  eigenfunctions  describe  the  natural  modes  of  vibration  of  the  fluid  and  the  eigenvalues 
determine  the  frequencies  at  which  these  mode  shapes  vibrate  (7:389). 

Now  consider  a  rectangular  domain  Q  such  that 

D  =  {x  =  (x,  y,  z)  :  0  <  x  <  IT,  0  <  y  <  L,  0  <  z  <  H}  .  (3.9) 

This  domain  is  a  specific  case  of  the  general  domain  discussed  above,  and  so  the  eigenfunc¬ 
tions  for  this  domain  can  be  found.  For  a  rectangular  domain,  the  spatial  operator  in  the 
wave  equation  is  the  Laplacian  operator  in  three-dimensional  Cartesian  coordinates,  so  for 
this  specific  domain,  the  homogeneous  wave  equation  is 

A  U  =  UXX  +  Uyy  +  UZZ  =  Utt  ■  (3.10) 

The  following  sections  go  through  the  process  of  using  this  eigenfunction  approach  by  first 
finding  the  time  component  in  the  same  way  it  was  found  above  and  then  by  finding  the 
eigenfunctions  of  the  Laplacian  operator. 

3.3  Dimensional  Analysis 

Dimensional  analysis  is  the  process  of  examining  the  dimensions  of  a  physical  quan¬ 
tity  to  describe  the  units  of  measurement.  In  Chapter  II,  it  was  stated  that  the  pressure 
is  simpler  to  measure  than  temperature  and  density  and  that  is  also  provides  more  accu¬ 
rate  results.  Furthermore,  a  microphone  measures  the  change  in  pressure  created  by  a 
sound  wave,  so  it  seems  logical  to  use  this  measurement  in  determining  the  effects  of  the 
environment  on  a  signal. 

Let  [a]  denote  the  dimension,  or  units,  of  the  quantity  a.  Let  u  (. x ,  t )  denote  the 
pressure  measured  at  location  x  at  time  t  and  define 

r  .  . .  F  KD  K 

[u(x,t)]  =  P=^  =  T^  =  T^  (3.11) 

where  K  is  the  unit  for  mass,  T  is  the  unit  for  time,  and  D  is  the  unit  for  length. 
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Since  the  differential  operator  is  a  rate  of  change  of  a  quantity  with  respect  to  another 
quantity,  the  dimension  of  the  derivative  is  the  original  dimension  divided  by  the  dimension 
of  the  quantity  the  derivative  is  taken  with  respect  to.  In  the  case  of  the  Laplacian 
differential  operator,  the  derivative  is  taken  with  respect  to  spatial  variables,  and  so 

[A«(x,f)]  =  -^  =  ^3  •  (3.12) 

Similarly,  the  dimension  of  the  quantity  \  utt  is 


1  V 

-o«tt  {x,t) 


1 


(T2\  f  P\  P  K 

[utt(x,t)\  —  I  Jy^  J  (  7^2  J  —  ~jyi  ~  rp2£)3  ' 


It  is  important  to  note  that  the  terms  A u  and  \ uu  have  the  same  dimension, 
physical  quantities  added  or  subtracted  together  must  have  the  same  dimension, 
this  in  mind,  it  is  easily  seen  that 


(3.13) 

Any  two 
Keeping 


f-Vu 


[Vu] 


K 

T2D 2 


implies  that 


1 

D' 


which  is  consistent  with  the  fact  that  the  atmospheric  attenuation  occurs  as  the  wave  travels 
through  space.  As  discussed  earlier,  atmospheric  attenuation  is  very  small  compared  with 
practical  boundary  attenuation,  and  so  examination  of  the  dimensionless  quantity  WF, 
where  W  is  the  largest  dimension  of  the  domain,  can  lead  to  a  justifiable  approximation 
of  r  =  0.  In  a  room  with  absolute  temperature  293  K  with  30%  relative  humidity,  the 
atmospheric  attenuation  is  only  .025  decibels  per  one  hundred  meters.  In  a  typical  room, 
this  attenuation  is  negligible,  and  so  the  dampening  term  discussed  in  Chapter  II  will 
be  dropped  from  the  wave  equation,  since  the  assumption  is  being  made  that  the  largest 
dimension  of  the  room  is  small  compared  to  the  speed  of  sound. 
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3-4  Duhamel’s  Principle 

3-4-1  Impetus.  When  solving  any  differential  equation,  ordinary  or  partial,  it  is 
always  easier  to  deal  with  a  homogeneous  problem  than  a  driven  one.  As  the  problem 
was  stated  earlier,  the  differential  equation  driven  by  a  forcing  function  with  its  associated 
initial  and  boundary  conditions  is 


Lu  = 

/ 

,  x  e  n, 

t  >  0 

u  ( x ,  0)  = 

0 

,  fed, 

c-+- 

II 

O 

ut  (x,  0)  = 

0 

,  x  €  fl, 

c-+- 

II 

O 

du  _ 

dn 

/ 3  (x)  u 

,  x  €  dfl, 

t  >  0 

while  what  is  sought  is  a  differential  equation  with  associated  initial  and  boundary  condi¬ 
tions  in  the  form 


Lu  = 
u  ( x ,  0)  = 

ut  (x,  0)  = 


0  ,  x  e  Q,  t  >  0 

0  ,  x  £  Q,  t  =  0 

g  ,  x  £  Q,  t  =  0 


(3.15) 


du 

dn 


(3  ( x )  u  ,  x  £  dll,  t  >  0 


There  is  a  theorem  that  accomplishes  this  goal,  called  Duhamel’s  Principle  (11:298). 
This  theorem  allows  for  the  reformulation  of  the  original  problem  into  a  homogeneous  one. 
This  homogeneous  problem  has  the  same  linear  differential  operator,  but  without  quiescent 
initial  conditions. 


3-4-2  Restatement  of  Problem.  If  u(x,t)  satisfies  the  original  differential  equa¬ 
tion,  then  v  ( x ,  t ;  r)  satisfies  the  initial  boundary  value  problem 


where 


Av 
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(3.16) 


(3.17) 
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Since  this  equation  is  homogeneous  in  a  rectangular  domain,  we  can  use  the  method  of 
separation  of  variables  to  solve  for  v(x,t;r),  which  will  yield  the  solution  in  terms  of 
u  (. x ,  t). 

3.5  Deriving  and  Inverting  the  RTF 

3.5.1  Separation  of  Variables.  The  geometry  of  a  domain  and  its  associated 
boundary  is  probably  the  largest  driving  factor  in  determining  the  form  that  the  solution 
of  a  differential  equation  will  take.  Using  a  rectangular  domain  for  the  wave  equation 
results  in  a  solution  in  the  form  of  infinite  series  of  trigonometric  functions,  while  solving 
the  same  equation  in  a  circular  domain  introduces  Bessel  functions  into  the  solution.  It 
is  fortuitous  that  most  rooms  are  roughly  rectangular,  since  this  allows  the  use  of  this 
method  to  solve  the  equation. 

The  method  of  separation  of  variables  is  a  valid  method  of  solving  a  partial  differ¬ 
ential  equation  in  rectangular  domains.  A  rectangular  domain  is  a  domain  in  which  the 
boundaries  with  respect  to  all  variables  are  not  dependent  on  any  other  variables.  For 
example,  in  Cartesian  coordinates,  a  circular  domain  is  invalid  for  using  separation  of  vari¬ 
ables  to  derive  a  solution  since  the  boundary  is  of  the  form  y  =  ±\/a2  —  x2.  However, 
performing  a  coordinate  transformation  into  polar  coordinates  will  result  in  a  rectangular 
domain  in  terms  of  the  radius  and  angle,  upon  which  a  separation  of  variables  method  can 
be  used.  Since  the  domain  in  question  is  already  rectangular,  the  equation  is  ready  to  be 
solved.  Begin  by  assuming  that  v  (x,  t;r)  =  Q  ( x ,  y,  z)  T  (f;  r). 

3.5. 1.1  The  Time  Component.  Assuming  that  v  takes  the  form  given  above, 
it  is  easily  seen  that  vtt  =  Q  (x,  y,  z)  T  (t;  r),  and  so  the  equation  becomes 

A  Q  _  f  _ 

Q(x,y,z )  c2T  (t;  r) 

where  oq  is  some  real  constant.  Focusing  only  on  the  second  and  last  parts  of  this  equation 
leads  to  the  ordinary  differential  equation 

T  —  a±c2T  =  0  (3.19) 
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which  can  take  on  three  forms,  each  based  on  the  value  of  a\.  These  three  cases  are 
described  next. 

Case  1:  ai  =  Xf  >  0.  In  this  case,  the  solution  is  the  linear  combi¬ 
nation  of  hyperbolic  functions 

T  (f ;  t)  =  A  sinh  (Aict)  +  B  cosh  (Ai ct)  .  (3.20) 

If  A  +  B  is  nonzero,  then  as  t  approaches  infinity  this  function  is  unbounded,  and  so  would 
violate  the  law  of  conservation  of  energy  for  the  system,  since  the  original  signal  has  finite 
energy.  Case  1  therefore  does  not  yield  a  nontrivial  solution. 

Case  2:  a\  =  A^  =  0.  In  this  case,  the  solution  is 

T(t;r)  =  At  +  B  .  (3.21) 

By  applying  the  first  initial  condition, 

v  (x,  t;t)  =  Q  ( x )  T  (r;  r)  =  0 

and  assuming  that  Q  (. x )  is  nonzero,  then 

T(r;r)  =  At  —  B  =  0 
B  =  At 

T(t;r)  =  A(t-T)  . 

And  by  applying  the  second  initial  condition, 

vt  (x,  r;  r)  =  -s  (r)  5  (x  -  x0) 

the  implication 

s  (r)  oc  At, 
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comes  up,  and  since  A  is  assumed  to  be  a  constant,  this  implies  that  the  signal  is  linear, 
and  so  violates  the  finite-energy  law  as  well.  Therefore,  Case  2  also  leads  to  only  a  trivial 
solution. 


Case  3:  a\  =  —  X±  <  0.  In  the  case  of  the  last  possibility,  the  solution 

is  of  the  form 

T  (t;r)  =  A  sin  (Xict)  +  B  cos  (Xict)  .  (3.22) 

Applying  the  first  initial  condition  and  assuming  that  Q  (x)  ^  0  for  all  x  yields 


v(x,t;r) 

=  Q  (x)  T  (t;t)  =  0 

(3.23) 

T  (t;  t) 

=  A  sin  (Aicr)  +  B  cos  (Aicr)  =  0 

(3.24) 

B 

=  —  Atan(Aicr) 

(3.25) 

T(t;r) 

=  A  (sin  (Ai ct)  —  tan  (Aicr)  cos  (Aicr)) 

(3.26) 

T(t;r ) 

=  C  sin  (Aic  (t  —  r)) 

(3.27) 

where  C  =  Acos(Aicr).  Since  this  is  a  constant,  this  is  allowed.  Because  the  second 
initial  condition  is  not  quiescent,  it  must  be  applied  after  the  eigenfunctions  have  been 
found.  As  this  is  a  finite  function  and  is  not  identically  equal  to  zero,  then  the  case  is 
valid,  and  so 

T  (t;  t)  =  A  sin  (Aic  (t  —  r))  (3.28) 


3.5. 1.2  The  Z  Component.  Now,  let  Q  ( x ,  y,z)  =  P  (x,  y )  Z  ( z ).  Then 


P'XX  +  Pyy 

P  (x,  y) 


7" 

A  \2 
-  —  -  Ai  —  a2 


(3.29) 


where  a 2  is  some  real  constant.  Examining  the  Z  component  first  yields  the  ordinary 
differential  equation 

Z"  +  (A?  +  a2)  Z  =  0  (3.30) 

which  has  characteristic  values 

±y/-x l -a2  ■  (3.31) 
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Again,  there  are  three  cases  which  will  describe  what  value  02  takes. 


Case  1:  a 2  <  —  Af. 


Let  ± 


\J  —  X2  —  02  =  ±q2.  In  this  case,  the 


quantity  under  the  radical  sign  is  positive,  and  so 


Z  (z)  =  Aeqz  +  Beqz  .  (3.32) 

On  the  floor,  at  z  =  0,  the  Robin  boundary  condition  is 

~Z'( 0)  =  f3zZ  (0)  (3.33) 

0  =  Z'(0)  +  PzZ(0)  (3.34) 

0  =  Aq  +  Bq  +  (3ZA  +  /3ZB  (3.35) 

B  =  -A  .  (3.36) 

On  the  ceiling,  at  z  =  H,  the  boundary  condition  is 

Z'(H)  =  (3zZ  (H)  (3.37) 

0  =  Z'  (H)  -  PzZ  {H)  (3.38) 

0  =  AqeqH  +  BqeqH  -  f3z(AeqH  +  BeqH)  (3.39) 

0  =  A(q-Pz)  +  B(q-Pz)  (3.40) 


which  leads  to  a  trivial  solution  since  A  =  —B.  Since  this  solution  would  lead  to  the 
entire  function  v  ( x ,  t;  r)  being  identically  zero,  this  case  is  invalid. 

Case  2:  02  =  —A^.  In  this  case  the  characteristic  root  is  simply,  in 
which  case  the  solution  takes  the  form 

Z  (z)  =  Az  +  B  ,  (3.41) 

and  so 

Z'  (z)  =  A  (3.42) 
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The  boundary  conditions  leads  to  the  following  two  equations: 


Z'(0)  +  /3*Z(0)  =  0 


(3.43) 


Z'(H)~PzZ(H)  =  0 


(3.44) 


but  it  is  clear  that  this  case  is  trivial  as  well,  and  so  case  2  is  also  invalid. 


values 


Case  3:  a2  >  —A2. 


In  this  case,  the  equation  has  the  characteristic 


±i\Ja,2  +  Af  , 


(3.45) 


the  solution  takes  the  form 


Z  (z)  =  A  sin  (qz)  +  B  cos  (qz) 


(3.46) 


where  q  =  \J <12  +  Af .  Using  the  boundary  conditions,  the  coefficients  can  be  solved  for: 


Z'(0)+/3*Z(0)  =  0  =  Aq  +  Bf3z 

q 


(3.47) 


B  =  -A- 


Pz 


Z'  (H)  —  P zZ  (H)  =  0  =  A  (q  cos  (qH)  —  [3Z  sin  ( qH ))  —  B  (q  sin  (qH)  +  (3Z  cos  ( qH )) 


tan  (qH)  = 


0  =  Asin(qH) 

2qPz 


+  A^—  sin  (qH)  +  Aq  cos  (qH)  (3.48) 


P 


q2- Pi 


(3.49) 


So  the  eigenvalues  associated  with  the  Z-cornponent  of  the  equation,  call  them  q  =  \i,  are 
the  intersections  of  the  curves  tan  (qH)  and  ,  and  so 

Q  P  z 


a2  =  —  A2 


(3.50) 
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The  eigenvalue-eigenfunction  solution  for  the  Z  component  is  therefore 


Z(z)  =  Al 


^sin  (Xiz) 


A; 

Pz 


cos  (A iz] 


1  =  1,2,3, ... 


3.5. 1.3  The  Y  Component.  One  the  left-hand  side  of  the  equation  remains 


the  expression 


Pxx  +  Pyy 
P(x,y) 


(3.51) 


Letting  P  ( x ,  y)  =  X  {x)  Y  (y)  results  in  the  elimination  of  the  last  partial  differential 
operators  and  leaves  only  the  ordinary  differential  equations 


X"  _  Y" 

A'(.r)  _V'(y) 

where  03  is  some  real  constant.  Looking  at  the  Y  (y)  portion  of  this  equation  leads  to  the 
ordinary  differential  equation 


\  ~  a3 


(3.52) 


Y"  +  (o3  +  A,2)  Y  =  0  . 


(3.53) 


Since  the  basic  ordinary  linear  differential  operator  is  acting  upon  Y  as  is  being  acted  upon 
Z,  the  same  arguments  for  the  value  of  the  constant  (03  +  Xf)  are  made  as  in  the  previous 
section.  So  by  using  the  same  methodology  used  for  the  Z  component,  the  Y  component 
solution  is 


Ym  (■ y )  =  Brn  ^sin  (A my)  -  ^  cos  (Xmy)^j  ,  m  =  1,2,3, ... 
where  Am  are  the  intersections  of  the  curves  tan  (qL)  and  ,  and 

1  Py 


(3.54) 


Q  —  A  m 


V a3  +  A;  . 


(3.55) 


3.5. 1-4  The  X  Component. 
component.  It  is  of  the  form 


X” 


The  final  component  of  the  equation  is  the  X 
a3X  =  0  (3.56) 
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Since  A§  is  a  positive  number,  the  equation  has  characteristic  roots 


My/a 3 


(3.57) 


with  the  associated  solution 

Xn  ( x )  =  C'n  ^sin  (Xnx)  -  ^  cos  (Anx)^  ,  n  =  1, 2, 3, ...  (3.58) 

where  An  =  —  y/ojj.  Using  all  of  the  relationships  between  the  Aj’s,  the  expression 

^1  =  =  \2  +  •  (3.59) 


3.5. 1.5  Joining  the  Components.  Now  that  the  four  components  of  the 
solution  have  been  found,  they  are  multiplied  together  to  solve  for  vij1Ujn(x,t-,T),  after 
which  the  principle  of  superposition  is  used  to  find  v  (x,  f;  r).  These  four  components  are 


,n  {t;  t) 

—  Al,m,n  sin  ( Xpm,n 

xn  (x) 

=  Cn  | 

^sin  (Xnx)  - 

Ym  (: V ) 

=  Dm 

^sin  (A my)  - 

Zi(z) 

=  Et  ( 

sin  (A iz)  -  - 

l,  77i,  n 

=  1,2, 

3,... 

A  n 

P~xcc 
A  m 
Py 
Xi 


(3.60) 

(3.61) 

(3.62) 

(3.63) 


Let  p  £  {x,y,  z},  D  £  {W,L,H},  and  o  £  {l,m,  n}.  Then  the  spatial  components 
of  the  solution,  the  eigenfunctions  with  their  associated  eigenvalues  for  each  component, 
can  be  generalized  by  the  function 


Fa(p;D) 


^sin  (A  ap) 


K 

Pp 


COS  (A ap) 


(3.64) 
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Factoring  out  the  constants  and  multiplying  the  generalized  functions  together  with  the 
time  component,  the  general  solution  is 

Vi,m,n  (x,  t;  t)  =  Ki,m,nFi  (z;  H)  Fm  (y;  L)  Fn  (x;  W)  sin  (\m,nc  ( t  -  r)) , 
v  (x,  t]T )  =  Kl,m,nFi  (z;  H)  Fm  (y; L)  Fn  (x;  W)  sin  (A z,m,nc  (t  -  r)) 

l,m,n 

where  Ki^m  n  =  A.i,m,nEiDrnCn. 

3.5.2  Determining  Coefficients.  Now  the  second  initial  condition  must  be  applied 
to  the  solution  to  determine  the  values  of  the  coefficients  K^m  n.  Recall  that  the  second 
initial  condition  is 


*(3.65) 
t  >{3-.j66) 


vt  (x,  t;  t)  =  -s  (r)  5  (x  —  xo)  (3.67) 

and  since  there  is  only  one  term  in  v  (x,  t;  r)  which  is  a  function  of  t,  the  derivative  is  easy 
to  compute: 


vt  (x,  tj  t)  —  ^  '  Ki,ni,nFi  (z;  .77)  F’yyj  (yj  _L)  Fn  (x;  VF)  Xi}rnnc cos  {Ximnc  (t  t))3.68) 

l,m,n 

vt  (x,  r;  r)  =  y  I<i,m,?iF  (z;  77)  Fm  (y;  L)  Fn  (x;  W)  Az,mjric  =  -s  (r)  (5  (x  -  x0X3.69) 

l,m,n 

Now,  by  taking  the  weighted  inner  product  over  the  domain,  the  coefficients  K]mn 
can  be  solved  for  algebraically.  Since  these  functions  are  all  linearly  independent,  then 
the  inner  product  has  no  weighting  function  by  the  Sturm-Liouville  Theorem  (7).  Then 
the  left-hand  side  of  the  inner  product  equation  is 

y  Ki,m,nFi  (z;  77)  Fm  (y;  L )  Fn  (x;  W)  Xi}m,nc  ,  F}/  (z;  77)  Fm/  (y;  L)  Fn>  (x;  IF) 

l,m,n 

(3.70) 

where  (it,  u)  =  J(1  u  (x)  v  (x)  dx.  Now,  by  taking  the  inner  product  over  the  domain,  the 
inner  product  is 
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fD 

/  Fa(p;D)  Fa,  (p;D)dp,  (3.71) 

Jo 

the  value  of  which,  since  the  eigenvalues  are  the  intersections  of  two  curves,  depend  on  the 
value  of  a  and  a'.  The  important  thing  is  that  the  eigenfunctions  are  orthogonal,  and  so 

for  a  ^  o’ ,  the  value  of  the  inner  product  is  equal  to  zero,  and  so  the  triple  sum  on  the 

left-hand  side  collapses  to  one  single  term 

(z,  l ',  H)  f  (y,  m!,  L)  f  (x,  n',W)  (3.72) 

/  (p,  a,  D )  =  [D  Fa  {p-  D)  Fa,  (p-  D)  dp  .  (3.73) 

Jo 

Now,  the  right-hand  side  of  the  inner  product  equation  is 

{-s  (r)  <5  (x  -  x0)  ,  Fv  (z;  H)  Fm,  (y;  L)  Fn,  (x;  W))  (3.74) 


and  since  the  Dirac  delta  function  acts  as  an  evaluation  when  it  is  in  the  integrand,  the 
right-hand  side  simplifies  to 

~8  (t)  (*o)  ,  (3-75) 

where 

V’ I,m,n  (£)  =  Fl  {z;  H)  Fm  (y;  L)  Fn  (x;  W)  (3.76) 


and  defining  fi.m:n  as 


fl,m,n=  f  (x,  n,  W)  f  (y,  m,  L)  f  (z,  l ,  H) 
the  equation  of  inner  products  is  now 


(3.77) 


c^l',m',n’Kl’,m’,n'fl’,m',n'  8  (t)  Ipi'  m'.n'  (*^o) 


(3.78) 


and  so, 


K, 


l,m,n 


~8  (r)  (^O) 

CA; 


(3.79) 
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3.5.3  The  Solution. 


3.5.3. 1  The  Time  Domain.  Now  that  the  initial  conditions  have  been 
utilized  to  solve  for  the  undetermined  coefficients,  the  solution  for  v  (x,  t;  r)  is 


v  (x 


r )  =  Y 

1,771,71 


^ l,m,n  (g)  j> l,m,n  <To)  H  (t)  shl  (A*,m,nC  (t  -  r))] 
c\i 

,m,nfl,m,n 


(3.80) 


and  so,  by  Duhamel’s  Principle,  the  time  domain  solution  is 

rt 


u 


(x,t)  =  y 

1,771,71 


^l,m,n  (g)  l,m,n  (fl))  fp  S  (t)  sin  (A;,m,nC  (t  ~  t))  dr 
cXl,m,nfl,m,n 


(3.81) 


3. 5. 3. 2  The  Laplace  Transform.  Now,  the  orginal  signal  s(t)  can  be  com¬ 
puted,  and  since  it  is  contained  within  an  integral,  a  transformation  must  be  applied  to 
eliminate  this  integral.  The  usual  transformation  of  choice  is  the  Fourier  transform,  but 
Fourier  transforms  deal  with  doubly  infinite  convolutions  (the  limits  of  integration  are 
positive  and  negative  infinity),  whereas  this  integral  has  finite  limits  of  integration.  The 
Laplace  transform  can  eliminate  this  integral  most  appropriately  (10:31),  and  so  taking 
the  Laplace  transform  of  both  sides  with  respect  to  the  time  variable  yields 


£ 


/  -.(r)sm(W(t-r))c*r}  =  -£  {s  (i)}  £  {sin  (Az)TO,nc  (t))}  (3.82) 

/»oo 

£{/(*)}  =  /  e-^f(t)dt  (3.83) 

Jo 


and  so  the  Laplace  transform  of  a  finite  convolution  is  the  product  of  the  Laplace  transforms 
of  the  elements  of  the  convolution.  Coupled  with  the  fact  that  £  {sin  (at)}  =  2+q2  this 
result  yields 


u  (x,  n) 
u  (x,  n) 


{h)  Y 

l, 771,71 

£  {u  (x,  t)} 


'4>l,m,n  O*')  ^ Pl,m,n  (To) 
{t2  +  Xj  ,7Tl,7l) 


(3.84) 

(3.85) 


Ht)  = 


(3.86) 
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3.5.4  The  Inversion. 


3.5.4- 1  The  Laplace  Transform.  The  Laplace  transform  of  the  original 
signal,  s  (/1)  can  now  be  solved  for  algebraically.  This  simple  operation  results  in  the 
Laplace  domain  solution  of  Problem  1: 


s{n) 


-u  (x,  n) 


E 


l’m’n  /l,m, „(#*»  + A?,m,n) 


(3.87) 


3. 5. 4-2  Back  To  the  Time  Domain.  Now,  by  taking  the  inverse  Laplace 
transform  of  both  sides,  the  time-dependent  and  source  location  independent  solution  to 
Problem  1  is 


s(t)  =  C  1 


—u  (xr,  n) 


5 


V,Z,m,nOrr)bi,m,„(To) 

fl,m,n  (/i2+A^mTl) 


> 


(3.88) 


where  xr  is  the  receiver  location.  Now,  this  Laplace  transform  will  most  likely  need  to  be 
solved  numerically  since  it  is  the  Laplace  transform  of  a  term  with  a  triple  infinite  sum  in 
the  denominator,  but  once  completed,  the  solution  to  Problem  1  is  complete. 


3.6  Filtering  Multiple  Sources 

3.6.1  Introduction.  In  this  next  section  Problem  2,  the  problem  of  filtering 
multiple  signals,  is  explored.  The  goal  of  this  section  is  to  generalize  Problem  1  to  filter 
extraneous  signals  from  a  recording.  In  this  section,  all  computations  are  made  in  the 
Laplace  transformed  time  domain,  and  so  all  results  will  be  stated  in  terms  of  the  time- 
dependent  signal’s  Laplace  transform,  unless  the  specified  equation  explicitly  contains  the 
variable  t.  In  this  section,  xr>{  denotes  the  location  of  the  ith  receiver  while  xs,i  denotes 
the  location  of  the  Ith  source. 

3.6.2  Restating  the  Problem.  Since  this  problem  concerns  multiple  signals,  the 
right-hand  side  of  the  modified  wave  equation  must  reflect  this  change.  Whereas  in 
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Problem  1  the  right-hand  side  is 


s  (' t )  5  (x-xq)  , 


(3.89) 


the  right-hand  side  of  the  wave  equation  in  Problem  2  has  a  sum  of  signals  at  distinct 
locations.  The  appropriate  expression  for  n  sources  at  location  i,  for  i  =  1,2,  ...,n  is 


Si  (*) 5  c 


x  —  X 


s,i ) 


(3.90) 


i=l 


If  there  are  also  multiple  receivers,  then  this  system  consists  of  the  convolution  equations 

.x  (*rj)E”=l^i,  m,n(xs,i)  fo  Si  (t)  sin  (Aj,  m,nc{t  ~  t))  dr 

u(xrJ,t)=  2^  - — - 7 -  (3.91) 

l.m.n 


i  =  1,2, ...,  n 


cAz 

j  =  1,2,  ...,m 


and  taking  the  Laplace  transform  of  both  sides  results  in  the  system  of  equations 

( xr,j )  X^=l  '4’l,m,n  ixs,i)  $i  {l1 ) 


u(xr,j,n)  =  ~Y^ 

l,m,n 


fl,m,n  (li2  +  \,m,n) 

u  (xrj,  fi)  =  C{u  (xrj,t)}  . 


(3.92) 

(3.93) 


In  a  more  compact  form,  the  system  can  be  written  as 


Uj(n)  =  -  y (n)  H  (xrj,xSji,  n) 


2—1 


H  {xr,j ,  xs,ii  fi)  —  Hij  —  2  ' 


'lPl,m,n  ixr,j)  '4>l,m,n  ( xs,i ) 


l  m  n  f l,tn,n  (li2  +  \,m,n ) 

Uj(n)  =  u{xr,j,n)  ■ 


(3.94) 

(3.95) 

(3.96) 


3.6.3  A  General  Solution.  Since  this  system  is  linear  with  two  separate  indices, 
it  can  be  written  as  a  matrix  equation,  where  m  is  the  number  of  receivers  and  n  is  the 
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number  of  sources.  The  m  x  n  matrix  of  transfer  functions  is 

Hu  H12  ■■■  Hln 

H-2i  H22  •  •  ■  H2n 

Hml  Hm2  •  •  •  Hmn 

and  the  system  of  linear  equations  is  the  matrix  equation 

His  =  u  (3.98) 

where 

si  (m) 

S2  (m) 

Sn  (/^) 

«1  (jj) 

U2  (jl) 

Um  (/i) 

This  matrix  equation  can  be  solved  to  find  the  vector  s.  As  is  the  case  with  all  linear 
equations,  certain  conditions  upon  BI  must  be  satisfied  in  order  for  there  to  be  a  solution, 
if  any.  BI  must  be  invertible,  which  is  discussed  in  a  later  section,  to  obtain  a  unique 
solution.  If  there  is  not  a  unique  solution,  as  in  the  case  where  the  number  of  sources 
outnumbers  the  number  of  receivers,  then  a  solution  of  minimum  norm  can  be  found,  and 
if  there  are  more  receivers  than  signals,  a  least-squares  solution  can  be  computed. 

3.6.4  A  Specific  Case-2  Receivers,  2  Sources.  To  apply  this  problem  to  a  real- 
world  problem,  consider  the  case  of  two  sources  and  two  receivers.  By  placing  no  limita¬ 
tions  on  the  location  of  the  receivers  or  sources,  this  example  confirms  that  a  single  signal 


(3.97) 
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can  be  found  by  the  received  signal  when  the  other  signal  is  considered  to  be  Sj  (//)  =  0. 
Since  this  is  a  linear  problem,  the  system  of  equations  reverts  to  the  original  problem,  but 
with  linearly  dependent  equations.  If  the  source  and  receiver  locations  are  constrained 
such  that  no  two  sources  or  receivers  are  collocated,  then  the  matrix  equation  can  be 
used  to  solve  for  the  individual  original  signals.  The  assumptions  from  Problem  1  are  still 
needed,  i.e.  the  sources  are  point  sources  not  located  on  the  boundary  and  the  microphones 
are  omnidirectional.  The  system  of  equations  for  this  problem  is 


Hu 

H\2 

si  (m) 

ill  (//) 

H21 

H‘2  ‘2  _ 

_  S2  (M)  _ 

_  U2  (n)  _ 

(3.101) 


Note  that  the  spatial  dependencies  of  s*  (jT)  and  Uj  (fa)  are  contained  in  the  subscripts. 
Now,  before  a  solution  to  this  problem  can  be  found,  the  nature  of  the  matrix  H  must  be 
determined,  as  it  determines  the  existence  of  a  solution. 


3.6.5  Limitations.  When  dealing  with  the  solvability  of  a  matrix  equation,  there 
are  three  possibilities.  The  first  is  that  there  exist  an  infinite  number  of  solutions.  In 
this  case,  the  number  of  sources  outnumbers  the  number  of  receivers  and  so  the  system  is 
under  determined.  When  this  is  the  case,  a  solution  which  has  minimum  norm  must  be 
sought,  most  easily  by  using  the  Moore- Penrose  pseudoinverse  of  the  form  where  m  <  n 
(4:114).  The  second  possibility  is  that  there  exists  no  solution  to  the  equation.  If  this 
is  the  case,  then  a  least-squares  solution  can  be  attempted  by  using  the  Moore-Penrose 
pseudoinverse  for  m  >  n.  The  final  possibility  is  that  there  exists  a  unique  solution  to 
the  system  of  equations.  In  this  case,  the  matrix  El  is  invertible,  and  standard  matrix 
algebra  can  be  used  to  solve  for  the  unknown  vector  s.  El  is  invertible  only  if  El  is  a  square 
matrix,  and  in  the  case  of  two  receivers  and  two  sources,  this  is  true,  but  the  existence  of 
a  solution  is  guaranteed  if  and  only  if  the  determinant  of  the  matrix  is  nonzero. 


3.6.5. 1  The  Determinant. 


Given  that  the  determinant  of  a  2  x  2  matrix  is 


a  b 
c  d 


ad  —  be  , 


(3.102) 
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the  determinant  of  the  matrix  HI  can  be  easily  computed  to  be 


det  (HI) 


H 11 

H 12 

h21 

h22 

=  huh22 


H12H21  . 


(3.103) 


Multiplying  any  arbitrary  two  transfer  functions  together  yields  the  sum 


rr  J r  ‘4’l,m,n  ( Xr,j )  'lPl,m,n  iXs,i )  V’p.q.r  ( Xr,k )  ixs,h) 

Uijtihk  ~  2^  7  /  2  T  \  2  7 

J  l,m,n  \H-  '  /'l,m,n ) 


l,m,n 


p,q,r 


fp,q,r  (fJ?  +  \,q,r) 


(3.104) 


and  with  a  bit  of  algebraic  manipulation,  the  expression  for  the  determinant  can  be  sim¬ 
plified  to 


M=  E 

l,m,n,pq,r 

Clearly,  if 


iXs,l)  iXs,2)  ^Pltm,n  (xr,l  )  V^q.r  (Xr, 2)  (^,2)  (®r,l)_ 

fl,m,nfp,q,r  (m  T  n  )  (m2  +  ^p,q,r) 

(3.105) 


^A/,m,n  (®r,l)  ^Pp,q,r  ( ®r,2 )  ^l,m,n  (xr,2)  V’p^r  (®r,l) 


(3.106) 


then  the  determinant  of  the  matrix  is  equal  to  zero.  This  will  happen  only  if  xr.i  =  xri2 
which  implies  that  if  the  receivers  are  located  in  the  same  position,  then  the  matrix  is  not 
invertible.  This  makes  sense  considering  that  if  two  microphones  are  collocated,  they  will 
record  the  same  signal.  There  would  be  no  difference  between  the  transformed  signals 
u  1  (/j)  and  u2  (//),  and  so  the  matrix  would  have  linearly  dependent  rows.  Similarly,  if  the 
sources  are  located  at  the  same  position,  then  the  determinant  would  be 


H  (xr!i  ,XSji,  n)  H  (xr}2,  XSji,  l-l)  -  H  (Xr72,  Xs,i,  /J,)  H  (xrti,XSji,  /i)  =  0 


(3.107) 


Since  H  (xrj,  xSji,  n)  has  no  dependence  upon  the  source  signal  (only  its  location), 
two  different  signals  propagated  from  the  same  location  will  be  impossible  to  separate  using 
this  method.  Therefore,  for  the  matrix  HI  to  be  invertible,  the  only  requirement  is  that 
the  two  sources  be  located  separately  and  the  two  receivers  be  located  separately.  If  a 
source  and  receiver  were  located  together,  then  there  is  no  problem.  In  fact,  this  would 
make  it  very  easy  to  filter  the  two  signals. 
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3. 6. 5. 2  The  Condition  Number.  The  determinant  is  not  the  only  factor  in 
the  solvability  of  a  system  of  linear  equations.  There  is  also  the  condition  number  of  the 
matrix.  The  condition  number  is  a  measure  of  the  stability  of  the  matrix  operator  in  terms 
of  solving  a  system  of  linear  equations.  If  a  matrix  is  nearly  linearly  dependent  (meaning 
that  the  determinant  is  close  to  zero),  then  the  matrix  will  have  a  high  condition  number. 
The  condition  number  is  not  unique;  it  depends  upon  which  norm  is  chosen  to  evaluate  the 
matrix  under.  Therefore,  there  are  a  number  of  different  possibilities  of  norms  that  can 
be  used.  When  dealing  with  matrices,  the  Frobenius  norm  is  a  good  norm  to  use  since  the 
evaluation  of  the  norm  requires  only  the  squaring  and  summing  of  the  matrix  components. 
The  Frobenius  norm  of  a  matrix  is  defined  as 


PH 2F  =  J2aij  •  (3-ios) 

Once  the  norm  is  chosen,  the  condition  number  of  the  matrix  A  is  defined  as 

Ka(A)  =  WA\\a\\A-1\\a  .  (3.109) 


where  a  is  the  choice  of  norm  used  in  finding  the  condition  number.  The  expression  can 
be  read  as  "The  condition  number  k  of  the  matrix  A  with  respect  to  the  norm  a  is  equal 
to  the  a  norm  of  A  times  the  a  norm  of  A-1."  As  stated  before,  the  Frobenius  norm  is 
useful  here  because  it  requires  minimal  calculations  and  the  inverse  of  a  2  x  2  matrix  has 
a  very  simple  formula.  The  Frobenius  norm  of  the  matrix  El  with  two  receivers  and  two 
sources  is  given  by 

PI|F  =  (i Ih+H212  +  H22l  +  H222y  .  (3.110) 

and  the  norm  of  HI-1,  where 


H11H22  —  H21H12 


H22  —H\2 

—H21  Hu 


(3.111) 
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is  given  by 


IT— 1  I 


1 

,-w  nun 


(H2 


11 


+  +  Hoi  +  Hoo) 2  — 


So  the  condition  number  of  the  matrix  El  is 


(3.112) 


kf  (El) 


PH  F 
det  (EQ 


(3.113) 


and  it  is  clear  that  the  condition  number  of  the  matrix  is  inversely  related  to  the  deter¬ 
minant  of  the  matrix  of  transfer  functions.  Simply  put,  if  the  receivers  or  sources  are 
located  closely  together,  then  the  condition  number  will  be  large,  and  the  matrix  will  be 
ill-conditioned,  yielding  a  potentially  unstable  result.  There  are  ways  around  using  a 
matrix  with  a  large  condition  number,  such  as  using  a  singular  value  decomposition  of  the 
matrix,  or  a  QR  decomposition,  but  these  measures  will  not  be  discussed  here. 


3. 6. 5. 3  The  Solution.  If  El  is  invertible  with  a  low  condition  number,  then 
the  system  of  equations  can  be  solved  simply  by  inverting  the  matrix  El.  Multiplying  both 
sides  of  the  equation  by  EM1  yields  the  solution 

H22Ul(u)  —  Hl2U2(u) 

det(H)  _  (3.114) 

HilU2{u)  —  H21Ul{T) 

det(H) 

3.1  Filtering  Sources  of  Variable  Location 

3.7.1  Restating  the  Problem.  In  Chapter  I,  Problem  3  was  introduced  as  the 
problem  of  filtering  the  signal  of  a  mobile  source.  In  Problems  1  and  2,  the  sources 
were  assumed  to  have  constant  location  during  the  recording,  and  so  the  position  was 
independent  of  time.  In  Problem  3,  however,  the  position  of  the  source  is  a  function  of 
time,  and  so  the  problem  is  somewhat  more  complicated  than  the  previous  problems.  The 
partial  differential  equation  needed  to  solve  this  problem  is  very  close  to  that  of  Problem 
1,  with  one  slight  change: 

A u - KUtt  =  s  (t)  5  (x  -  x0  (t))  ,  (3.115) 

cz 


si  Qt> 
S2  (n) 
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where  xq  (t)  is  the  location  of  the  source  at  time  t.  Applying  Duhamel’s  Principle  to  this 
problem  yields  the  same  homogeneous  differential  equation,  but  with  the  initial  condition 

vt  (x,t;t)  =  -s  (r)  S  (x  -  xo  (r))  .  (3.116) 

The  analysis  of  the  problem  is  very  much  the  same  as  it  was  with  Problem  1,  but 
when  the  weighted  inner  product  is  taken  over  the  domain,  the  right-hand  side  of  the  inner 
product  equation  is  somewhat  different: 

{-s  (r)  5{x-x o  (r))  ,  Fv  (z;  H)  Fm,  (y,  L)  Fn,  (x;  W))  (3.117) 

Whereas  in  Problem  1,  only  the  time-dependent  signal  was  a  function  of  r,  any  term  that 
is  dependent  upon  the  source  location  xq  has  a  time-dependency.  The  right-hand  side  of 
the  inner  product  equation  is 


-s(r)V,i',m',n'(*o(r))  •  (3.118) 

The  complications  do  not  manifest  until  the  Laplace  transform  of  the  finite  convolu¬ 
tion  is  taken.  Since  the  source  location  is  dependent  upon  time,  it  cannot  be  factored  out 
of  the  integral,  and  so  the  Laplace  transform  is  more  complicated.  The  convolution  is 

(r)  V’ i',m',n '  (*o  (t))  sin  (A itm,nc  (t  -  r))  dr  .  (3.119) 

To  simplify  this  convolution,  define 

G  (t)  =  -S  (t)  (*0  ( t ))  (3.120) 


Then  the  integral  becomes 


t))(It  , 


(3.121) 
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and  by  taking  the  Laplace  transform,  the  expression  is  similar  to  that  of  Problems  1  and 
2,  but  with  G  (n)  instead  of  s  (/i): 


C\  I  G  (r )  sin  {\i,m,nc  (r  —  t))dr  f  =  G  (/i)  C  {sin  (3.122) 


G(»)=£{G(t)}  . 


(3.123) 


G  (/j,)  can  be  expressed  in  the  same  way  as  s  (/j)  in  the  previous  problems.  So, 


G(fj)  = 


u  (x,  n) 


,m,n  (z)V’i  ,m,n  (gp)  ’ 

1,171,11  /l,m,n(M2  +  (Al,m,„c)2j 


E 


(3.124) 


but  a  solution  for  s  (ji)  is  desired,  so  some  analysis  on  G  (//)  must  be  performed.  By 
defining 

M  (t)  =  {xo  (t))  ,  (3.125) 

then  taking  the  Laplace  transform  of  both  sides  of 


G  (t)  =  s  (' t )  M  (■ t ) 


(3.126) 


yields 

POO 

G  (n)  =  s  (n)  *  M  (//)  =  /  s  (jj)  M  {n  —  rj)  drj  (3.127) 

4o 

since  the  Laplace  transform  of  a  product  of  two  functions  is  the  convolution  of  their  Laplace 
transforms.  By  taking  a  second  Laplace  transform  with  respect  to  the  transform  variable 
r),  the  convolution  is  converted  to  a  multiplication  of  the  two  doubly  transformed  functions, 
and  s  (t)  can  be  solved  for  by  taking  the  inverse  Laplace  transforms  of  both  sides,  first  with 
respect  to  rj,  and  second  with  respect  to  /_/. 

There  is  a  simpler  way  to  achieve  an  approximation  of  this  result,  however.  Assuming 
that  the  velocity  of  the  source  is  small  compared  to  the  speed  of  sound,  that  is 

<  c,  (3.128) 


dx 

dt 
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then  by  sampling  the  source  location  over  small  time  intervals  and  summing  the  signals 
over  the  length  of  the  recording,  the  problem  can  be  converted  to  one  similar  to  Problem 
1,  but  with  the  signal  computed  on  each  interval  and  then  summed.  In  essence,  Problem 
1  is  applied  to  each  interval.  To  achieve  this  goal,  let 

xq  ( U )  =  bi,  ti- 1  <t<ti,  *  =  1,2, ...,  N  (3.129) 


where  to  =  0  and  t^  = 


m  (x, t )  = 

u  (. x ,  t)  = 


T,  the  total  length  of  time  of  the  recording. 


E 

l.m.n 


'I’l.m.n  (®)  (hi)  fo~  Si  (r)  sin  (Xt} 


Then,  for  ti- 1 
m,n  (t-r))  dr 


N 

^2 Ui  *) 

i=l 


(3.130) 

(3.131) 


Thus,  applying  the  results  of  Problem  1  to  each  interval  of  the  recording  yields  the  original 
signal  for  that  interval,  and  so  the  total  signal  is  simply  the  sum  of  the  signal  over  each 
interval: 

N 

«(*)  =  £*(*)  (3.1.32) 

i= 1 

For  these  definitions  of  u  (x,  t )  and  s  (t)  to  make  sense,  then  outside  of  the  interval 
then  Ui(x,t )  and  Si(t)  must  be  identically  equal  to  zero.  Continuity  of  the 
functions  and  their  derivatives  at  the  endpoints  of  the  integral  must  also  be  considered  so 
that  the  signal  is  smooth  after  being  filtered.  While  the  number  of  calculations  required  to 
perform  this  is  N  times  the  number  of  calculations  required  to  filter  a  source  with  constant 
location,  it  is  more  applicable  to  the  requirements  of  filtering  the  signals  from  a  room  of 
people.  Efficiency  must  be  considered  when  deciding  how  many  points  to  sample,  as  well 
as  the  difficulty  in  determining  the  source  location  at  any  given  time. 


3.8  Conclusion 

3.8.1  Results.  In  Problem  1,  it  is  shown  that  the  linear  differential  operator 
acting  on  a  function  u  (x,  t )  with  a  driving  term  s  (t)  5  (x  —  To)  is  invertible,  and  the  original 
signal  can  be  found  after  it  has  been  acted  on  by  the  environment  in  which  it  propagates, 
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as  long  as  the  geometry  of  the  domain  is  sufficiently  simple.  In  this  case,  the  solution 
was  found  for  a  rectangular  room,  but  it  was  shown  in  Section  2  that  the  problem  can 
be  generalized  easily,  where  the  only  difference  in  the  procedure  is  that  of  finding  the 
eigenfunctions.  Rectangular  domains  are  especially  good  for  solving  analytically  due  to 
the  fact  that  the  Laplacian  operator  for  this  coordinate  system  is  separable  in  the  spatial 
domains.  Because  of  this,  the  eigenfunctions  were  relatively  easy  to  find.  However,  with 
different  domains,  the  eigenvalues  and  eigenfunctions  can  be  harder  to  compute,  such  as 
the  case  with  a  cylindrical  domain.  The  Laplacian  operator  in  a  cylindrical  Cartesian 
domain  is  not  separable,  except  in  the  z-direction.  However,  after  applying  a  coordinate 
transformation  from  Cartesian  coordinates  to  polar  coordinates,  the  operator  is  separable 
with  the  difference  being  that  the  coefficients  are  no  longer  constant.  Because  of  this 
difference,  functions  containing  Bessel  functions  emerge  as  eigenfunctions.  Furthermore, 
in  some  irregular  domains,  a  closed-form  solution  for  the  eigenfunctions  may  not  even  be 
possible,  and  so  a  numerical  method  such  as  a  finite  element  method  would  be  required  to 
find  the  eigenfunctions.  So,  it  is  easy  to  see  how  a  change  in  the  domain  geometry  leads 
to  different  mode  shapes  for  the  room  transfer  function. 

In  Problem  2,  the  circumstances  under  which  two  signals  can  be  filtered  were  ex¬ 
plored,  and  the  results  were  applied  to  the  case  of  two  receivers  and  two  sources.  While 
this  is  the  simplest  case  conceivable,  it  is  not  very  practical.  It  is  easy  to  generalize  the 
solution  under  two  sources  and  receivers  to  a  larger  problem  simply  by  expanding  the 
dimensions  of  the  matrix  HI,  but  this  leads  to  more  computationally  intensive  equations. 
The  circumstances  under  which  the  HI  matrix  is  ill-conditioned  were  also  discussed,  as  well 
as  the  effects  on  the  solution  when  HI  is  ill-conditioned,  leading  to  a  potentially  unstable 
result.  This  instability  can  be  avoided  by  using  singular  value  decompositions  or  QR 
decompositions. 

In  Problem  3,  a  source  with  variable  location  was  discussed.  The  difficulty  in  solving 
this  equation  comes  about  from  the  extra  terms  left  in  the  convolution  integral,  and  the 
method  of  applying  a  second  Laplace  transform  was  offered  as  a  means  of  extracting  the 
original  signal.  A  sampling  method  was  then  devised  that  yields  an  approximation  of 
the  total  signal  by  sampling  the  source  location  during  a  time  interval  and  summing  the 


3-26 


solutions  on  these  intervals  together,  while  assuming  that  the  velocity  of  the  source  is 
small  with  respect  to  the  speed  of  sound.  While  more  costly  in  execution  than  Problem 
1,  it  is  far  more  practical  than  considering  a  human  source  stationary  during  the  time  of 
recording. 

3.8.2  Practicality.  As  discussed  before,  the  largest  obstacle  to  overcome  in  solving 
Problem  1  for  a  particular  domain  is  the  geometry  of  the  room.  A  different  domain  could 
produce  vastly  different  results.  The  possibility  of  not  knowing  the  parameters  of  the 
room  must  also  be  considered,  such  as  the  coefficients  of  refiection  and  whether  or  not  the 
boundaries  transmit  sound  outside  of  the  domain. 

The  question  of  real-time  filtering  is  important  since  the  importance  of  the  informa¬ 
tion  being  filtered  could  be  very  high.  In  the  case  of  intelligence  gathering,  decisions  may 
be  based  on  the  information  contained  in  the  recordings,  and  so  it  is  necessary  that  this 
information  be  found  as  quickly  as  possible,  so  real-time  filtering  is  somewhat  important 
here.  In  the  case  of  using  a  filter  to  eliminate  noise  heard  in  a  phone  call  from  a  loud 
room,  real-time  filtering  is  a  necessity  for  the  product  to  be  marketable  to  the  public. 
However,  in  the  case  of  recording  boardroom  meetings,  the  need  for  real-time  filtering  is 
nonexistent,  since  the  filtered  signals  would  be  stored  for  future  use  anyway.  There  is  no 
inherent  problem  in  the  mathematics  of  the  idea  of  real-time  filtering  using  this  method, 
since  the  expression  for  u  (x,  t )  is  not  dependent  upon  the  total  duration  of  the  recording, 
T.  Therefore,  a  computer  program  could  filter  each  second  recorded  in  a  period  of  time 
dependent  only  upon  the  number  of  signals  and  number  of  receivers,  and  of  course  the 
speed  at  which  it  can  perform  calculations.  The  implementation  of  such  a  program  is 
developed  and  discussed  in  Chapter  IV. 

The  final  limitation  on  the  practicality  of  this  method  is  that  of  being  able  to  de¬ 
termine  source  locations.  It  is  entirely  possible  to  determine  source  locations  using  an 
echolocation  procedure  much  the  same  as  bats  use  ultrasonic  sonar  to  determine  the  loca¬ 
tion  of  their  prey  (source  here).  Whether  or  not  this  active  form  of  source  localization  or  a 
passive  form  using  the  difference  in  times  at  which  the  receivers  actually  receive  the  signal 
depends  on  the  circumstances  under  which  the  recording  is  being  made.  If  it  is  covert  or 
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would  interfere  with  the  conversation  being  recorded,  then  a  passive  form  of  source  local¬ 
ization  must  be  used.  If  it  is  a  meeting  that  everyone  knows  is  being  recorded  or  where 
there  is  no  way  that  they  could  know  about  the  active  localization,  then  an  active  form  of 
source  localization  can  be  used. 
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IV.  Application  Analysis 


4-1  Introduction 

This  chapter  attempts  to  determine  the  net  effect  of  the  environment  on  a  signal 
in  a  computational  approach  using  the  results  of  Chapter  III.  The  lack  of  experimental 
data  severely  hinders  the  full  extent  of  the  application  of  these  results,  but  the  procedure 
can  be  applied  to  known  signals  to  determine  how  well  the  assumptions  made  in  solving 
Problems  1  and  2  fit  with  the  real-world  effects  of  an  environment  on  acoustic  signals. 
Of  the  effects  measured  computationally,  the  net  change  in  signal  content  in  the  Laplace 
transform  domain  will  be  measured  on  single-frequency  tones.  Without  more  powerful 
computing  tools,  the  full  analytic  value  of  these  approximations  is  diminished,  but  it  serves 
to  show  the  number  of  calculations  required  to  filter  the  signal  accurately.  For  example, 
the  loop  in  the  program  that  determines  the  components  of  the  matrix  HI  can  require 
millions  upon  millions  of  cycles,  depending  on  how  accurate  the  results  are  needed  to 
be.  In  the  computations  performed  in  this  chapter,  the  loop  that  found  these  components 
cycles  approximately  two  million  times  to  achieve  a  numerical  approximation  of  the  transfer 
function. 

The  program  used  to  implement  the  solutions  found  in  Chapter  III  is  MATLAB, 
version  7.0.1.  The  original  signal  is  input  symbolically  and  then  quantified  numerically 
to  aid  in  computation  speed.  The  original  program  used  the  symbolic  algebra  toolbox  for 
MATLAB  almost  exclusively,  but  the  computing  power  needed  to  perform  this  symbolically 
is  very  high,  and  the  times  for  computing  the  transfer  functions  alone  in  this  manner  took 
more  than  15  minutes  per  calculation.  Since  this  method  was  so  inefficient,  it  was  dismissed 
as  a  possibility  for  finding  the  transfer  functions.  The  numerical  approach  greatly  reduces 
the  computational  time  needed,  but  still  takes  approximately  6  minutes  to  implement  for 
a  simple  single-frequency  tone  and  so  would  be  very  inefficient  for  computing  the  transfer 
function  for  a  highly  intricate  pattern  of  a  typical  human  voice  speaking  for  an  extended 
period  of  time.  The  actual  MATLAB  code  used  in  determining  these  transfer  functions 
is  in  Appendix  A. 
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Because  the  computational  power  needed  for  computing  a  single  transfer  function  is 
so  high,  the  computation  of  multiple  transfer  functions  is  avoided,  but  is  discussed  in  the 
Section  3  of  this  chapter.  The  MATLAB  code  to  implement  the  two  receiver/two  source 
case  is  also  contained  in  Appendix  A.  Suggestions  for  how  to  implement  this  code  more 
efficiently  are  contained  in  Chapter  V. 

Finally,  in  all  calculations,  the  atmospheric  attenuation  was  ignored  since  the  room  is 
assumed  to  be  sufficiently  small  to  ignore  the  0.025  rn  attenuation  due  to  atmospheric 
friction.  This  eliminated  some  unnecessary  complications  in  the  writing  of  the  code.  To 
eliminate  further  complications  in  the  code,  the  eigenvalues  were  approximated  by  for 

a  =  1,  2,  3, .  It  was  shown  in  Chapter  III  that  the  eigenvalues  were  the  intersection  of  a 

rational  function  and  a  tangent  function,  but  for  low  (3  and  low  index,  this  approximation 
is  fairly  valid. 

4-2  Single  Frequency  Tones 

The  human  auditory  range  in  young  adults  is  20  Hz  to  20  kHz.  This  section  uses  the 
MATLAB  code  for  one  source  and  one  receiver  to  determine  the  effect  of  an  environment 
on  a  single- frequency  tone.  The  parameters  for  the  environment  that  were  used  are 

W  =  10  m 

L  =  10  m 

H  =  3  m 

xr  =  (2  m,  2  m,  2  m) 

xs  =  (9  m,  9  m,  1.5  m) 

Px  =  io 

Py  =  10 

Pz  =  10  • 

The  room  is  modeled  after  a  typical  medium-sized  room  with  sufficiently  small  coefficients 
of  refiection  to  avoid  extraneous  signal  reflection  in  order  to  aid  in  computation. 
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4-2.1  Low  Frequency  Range-20  Hz.  In  the  very  lowest  range  of  human  auditory 
perception,  the  change  in  the  Laplace  transform  of  the  signal  is  greatest.  The  two  changes 
were  computed  at  two  distinct  source  locations-one  at  a  moderate  distance  from  the  re¬ 
ceiver  and  the  other  very  close  to  the  receiver.  The  lack  of  accuracy  in  these  computations 
is  intuitively  seen  in  the  fact  that  the  change  in  the  Laplace  transform  of  the  original  and 
received  signals  is  identical.  This  change  in  the  signal  is  computed  by  approximating  the 
integral 

r 20000  r-20000 

1=  {s  (/a)  -  u  (/a))2  d/a  =  s  {/a)  (1  -  H  (n))2  dn  (4.1) 

J-20000  .1—20000 

with  the  rectangular  rule  of  Riemann  integration.  The  width  of  each  rectangle  is  10,  and 
so  the  total  number  of  terms  in  the  Riemann  sum  is  4001. 

4-2. 1.1  Source  Location  1.  Using  the  source  location  given  above,  which 
places  the  source  at  a  large  distance  from  the  receiver  with  respect  to  the  size  of  the  room, 
the  original  signal  was  assumed  to  be 

s  ( t )  =  sin  (27 r  (20)  t) 

which  oscillates  fairly  rapidly-20  cycles  per  second.  The  Laplace  transform  of  this  signal 
is  shown  in  Figure  4.1. 

The  transfer  function  is  determined  solely  by  the  environment  in  which  the  signal 
propagates,  and  so  if  the  parameters  of  the  room  remain  constant  when  two  different  signals 
are  being  propagated  from  the  same  location,  then  the  two  signals  are  acted  upon  by  the 
same  transfer  function.  Since  we  are  simply  looking  at  the  transformed  transfer  function 
along  the  real  line,  and  since  there  are  poles  along  the  imaginary  axis  of  the  transform 
domain,  there  is  a  peak  in  the  value  of  the  transfer  function  close  to  n  =  0.  The  transfer 
function  for  the  room  described  above  is  shown  in  Figure  4.2.  The  scale  on  Figures  4.1  and 
4.2  are  very  different,  but  when  scaled,  both  are  very  nearly  zero  except  when  close  to  the 
origin.  The  received  signal  can  be  found  simply  by  multiplying  the  transformed  signal  and 
the  transformed  transfer  function  together.  The  received  signal  in  the  transform  domain 
is  shown  in  Figure  4.3. 
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Figure  4.1  Laplace  Transform,  20  Hz. 
Source  Location  2.  The  second  source  location  is 


xs  =  (-1  m,  9  m,  2.9  m) 


Since  the  signal  is  the  same,  only  the  transfer  function  changes  with  respect  to  the  transform 
variable  ^ i  from  the  previous  source  location.  The  transfer  function  for  this  source  location 
is  shown  in  Figure  4.4. 

4-2. 1.2  Middle  Frequency  Range-2  kHz.  The  middle  range  of  human  au¬ 
ditory  perception  is  around  2  kHz,  which  oscillates  much  more  rapidly  than  the  20  Hz 
signal,  and  so  to  model  it  accurately  in  the  time  domain  requires  the  use  of  approximately 
100  times  more  sample  points  than  with  the  20  Hz  signal.  Therefore,  it  is  more  useful  to 
perform  these  calculations  completely  in  the  Laplace  transform  domain. 


be 


Source  Location  1.  In  this  section,  the  original  signal  is  assumed  to 


s  (t)  =  sin  (27t  (2000)  t) 


(4.2) 
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Figure  4.2  Transfer  Function,  Laplace  Domain 


which  has  a  Laplace  transform  of 


~  /  40007T 

s  (N  =  — - o  ■ 

/j,2  +  (40007r)2 


(4.3) 


the  graph  of  which  is  contained  in  Figure  4.6. 

Whereas  the  20  Hz  signal  dropped  off  closer  to  zero,  this  signal  drops  off  much 
farther  from  the  origin.  Furthermore,  the  maximum  value  of  the  transform  is  on  the  order 
of  hundredths  whereas  the  maximum  value  of  the  20  Hz  signal  transform  was  on  the  order 
of  1.  The  transformed  received  signal  is  given  in  Figure  4.7. 


Source  Location  2.  In  this  location,  the  received  signal  from  the 
original  2  kHz  signal  is  much  the  same  as  that  for  source  location  1,  since  the  transfer 
functions  are  not  so  different  for  these  two  source  locations. 


4-  2. 1.3  High  Frequency  Range-20  kHz.  At  the  upper  frequency  range  of 
human  auditory  perception,  the  transformed  signal  looks  much  the  same  as  that  of  the  2 
kHz  signal,  but  with  more  curvature.  As  with  the  2  kHz  signal,  the  number  of  sample 
points  for  any  time-domain  computations  is  prohibitive,  and  so  only  the  Laplace  transform 
domain  will  be  dealt  with.  The  Laplace  transform  of  the  20  kHz  signal  is  given  in  Figure 
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Figure  4.3  Received  20  Hz  Signal,  Transform  Domain. 

4. 8. By  applying  the  transfer  function  to  this  transformed  signal,  the  received  signal  is 
attained,  which  is  shown  in  Figure  4.9. 

4-2.2  Multifrequency  Signals.  Different  multifrequency  signals  can  be  considered 
for  implementation  in  the  MATLAB  program,  but  the  received  signals  cannot  be  computed 
accurately  due  to  their  complexity.  There  are  two  different  types  of  multifrequency  signals 
that  can  be  implemented  easily  by  modifying  the  program,  chirps  and  simultaneous  tones 
of  different  frequency. 

4-2.2. 1  Chirps.  Chirps  are  signals  that  either  increase  or  decrease  in  fre¬ 
quency  with  time.  The  rate  of  increase  or  decrease  affects  the  transformed  signal  by  its 
Laplace  transform.  One  example  of  an  increasing  tone  that  ascends  in  pitch  quickly  is 
the  time-signal  s(t)  =  sin  (27rf2)  .MATLAB  encounters  difficulty  in  attempting  to  find  a 
closed-form  solution  for  the  Laplace  transform  of  a  chirp,  and  so  it  must  be  computed  nu¬ 
merically.  The  graph  of  the  numerically  computed  Laplace  transform  is  shown  in  Figure 
4.10.  It  is  apparent  that  there  is  a  pole  at  the  origin,  and  so  when  the  transfer  function  is 
applied  to  this  transformed  signal,  the  pole  will  be  carried  through  to  the  received  signal 
by  multiplication. 
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Figure  4.4  Transfer  Function,  Laplace  Domain,  Source  Location  2. 

An  example  of  a  slowly  ascending  chirp  is  easy  to  find.  Simply  decrease  the  value  of 
the  exponent  in  the  rapidly  ascending  chirp,  but  without  going  to  or  below  1.  By  lowering 
the  exponent  to  only  1.5,  the  chirp  ascends  in  pitch  much  more  slowly,  as  shown  in  Figure 
4.11. 

Since  the  frequency  increases  much  more  slowly  than  s  (t)  =  sin(27rt2),  even  with 
only  a  slight  decrease  in  the  exponent,  MATLAB  has  much  less  difficulty  in  computing 
the  Laplace  transform  analytically,  and  with  a  program  that  computes  the  received  signal 
well,  will  yield  a  better  result  than  with  a  rapid  ascension  chirp. 

Descending  chirps  behave  very  similar  in  MATLAB  to  their  corresponding  ascending 
chirps.  By  simply  replacing  the  exponent  of  the  time  variable  with  its  reciprocal,  the 
chirp  will  descend  in  pitch  at  the  same  rate  that  its  counterpart  ascends.  A  phase  shift 
must  be  introduced  as  well,  to  avoid  a  singularity  in  the  derivative  at  t  =  0.  The  rapidly 
ascending  chirp’s  descending  counterpart  is  s  (t)  =  sin  ( 2iry/t  +  l) ,  which  decreases  much 
more  slowly  than  its  counterpart  increases.  Note  that  the  scale  on  the  t- axis  of  Figure 
4.13  is  ten  times  that  of  Figure  4.11. 


4-7 


Figure  4.5  Laplace  Transform,  2  kHz  Signal. 


The  counterpart  to  the  slowly  ascending  chirp  is  the  descending  chirp 

s  (t)  =  sin  ^27t  \J(t  +  l)2^  (4.4) 

Note  again  that  the  descending  chirp  descends  much  more  slowly  in  pitch  than  its  counter¬ 
part,  as  Figure  4.14  has  a  scale  ten  times  that  of  Figure  4.10.  Like  the  rapidly  ascending 
chirp,  MATLAB  encounters  difficulty  in  performing  the  Laplace  transform  of  this  signal 
because  of  its  high  initial  frequencies. 

4-2. 2. 2  Simultaneous  Tones  of  Different  Frequency.  The  touch-tone  phone 
system  uses  a  signal  consisting  of  two  distinct  tones  overlapped  for  each  button  on  a 
telephone.  One  component  of  the  signals  for  the  buttons  in  a  row  all  have  the  same  tone 
and  one  component  of  the  signals  for  the  buttons  in  a  column  share  the  same  tone.  This 
is  an  example  of  a  use  for  a  ditonal  frequency.  A  human  voice  is  comprised  of  many 
tones  overlapping  and  being  propagated  simultaneously,  so  these  signals  are  the  closest 
to  human  voices  of  those  that  have  been  analyzed.  An  example  of  a  ditonal  signal,  i.e. 
a  signal  comprised  of  two  tones  is  s  ( t )  =  sin  (2-7T  (200t))  +  sin  (27 r  (400i)),  whose  graph  is 
displayed  in  Figure  4. 15. The  Laplace  transform  of  this  signal  is  much  easier  for  MATLAB 
to  compute  than  that  of  chirps  of  any  kind,  since  this  is  simply  the  sum  of  two  single-tone 
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Figure  4.6  Laplace  Transform,  2  kHz  Signal. 


frequencies,  and  the  Laplace  transform  is  a  linear  transformation.  In  fact,  given  the  value 
of  the  amplitudes  at  given  frequencies,  a  fair  approximation  of  a  human  voice  could  be 
made  using  this  program.  Once  again,  however,  the  filtering  of  such  a  signal  is  impossible 
with  the  program  provided  in  Appendix  A.  The  graphs  of  signals  with  multiple  tones, 
however,  becomes  much  more  complicated  than  that  in  Figure  4.15.  For  example,  a  signal 
consisting  of  four  frequencies,  200  Hz,  850  Hz,  1700  Hz,  and  2500  Hz  with  amplitudes  of  2, 
5,  3,  and  4,  respectively  (these  numbers  are  chosen  somewhat  randomly,  only  constrained 
to  fall  within  the  low  to  middle  range  of  human  auditory  perception)  has  the  signal 

4 

s(t)  =  y^Asin  (27 T(fjt))  (4.5) 

2=1 

Ai  =  (2, 5, 3, 4)  (4.6) 

fi  =  (200,  850,  1700,  2500)  (4.7) 

which  has  a  graph  displayed  in  Figure  4.15.  The  complexity  of  the  signal  is  very  clear,  and 
this  is  only  for  four  distinct  tones.  With  the  continuous  range  of  tones  contained  in  the 
average  human  voice,  the  signal  only  gets  much  more  complicated,  making  it  very  difficult 
to  filter  in  a  computer  program. 
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Figure  4.7  Laplace  Transform,  20  kHz. 


Figure  4.8  Transformed  Received  Signal,  20  kHz. 
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Figure  4.9  Rapid  Ascension  Chirp,  Time  Domain. 


Figure  4.10  Laplace  Transform,  Rapid  Ascension  Chirp. 
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Figure  4.15  Four- Tone  Signal  of  Varying  Amplitude. 


4-14 


V.  Conclusions  and  Recommendations 


5.1  Summary 

The  initial  problem  of  finding  and  inverting  a  room’s  acoustic  transfer  function  begins 
with  finding  the  rule  by  which  acoustic  signals  propagate  and  then  by  finding  the  rule 
by  which  the  propagation  affects  the  signal.  Once  this  rule  is  found,  the  effects  of  the 
environment  can  be  reversed,  thereby  yielding  the  original  signal.  When  multiple  receivers 
and  multiple  sources  are  involved,  the  set  of  effects  on  the  signals  is  described  by  a  system 
of  linear  equations  in  a  transformed  domain,  which  can  be  expressed  in  terms  of  a  matrix. 
Under  certain  conditions,  this  matrix  can  be  inverted,  and  a  unique  solution  for  the  set 
of  signals  can  be  found,  but  quite  often  these  conditions  are  not  met,  resulting  in  a  least- 
squares  or  minimum  norm  solution.  On  top  of  these  constraints  is  the  fact  that  the 
computations  are  made  with  the  assumption  that  the  sources  are  in  a  fixed,  known  location 
throughout  the  duration  of  the  recording.  Once  this  assumption  is  false,  the  problem  of 
inverting  even  a  single  transfer  function  becomes  much  more  difficult,  involving  multiple 
Laplace  transforms. 

With  all  of  these  conditions  met,  it  is  theoretically  possible  to  filter  acoustic  signals 
exactly,  but  the  implementation  of  this  process  in  a  program  as  advanced  as  MATLAB  can 
be  prohibitive,  even  for  simple  signals  such  as  single-frequency  tones.  When  considering  the 
acoustic  complexities  of  the  human  voice,  the  implementation  could  take  a  prohibitively 
long  period  of  time.  The  goal  of  a  real-time  filter  is  one  that  would  benefit  both  the 
Department  of  Defense  and  private  industry.  With  real-time  filtering,  information  from  a 
bugged  room  could  be  processed  immediately,  yielding  potentially  life-saving  information, 
and  conversations  held  over  a  mobile  phone  in  a  noisy  room  would  not  be  as  frustrating  as 
they  usually  are.  Products  such  as  the  adaptive  earpiece  Jawbone  take  advantage  of  this 
technology  to  filter  the  noise  from  a  loud  room  (1). 

The  raw  computing  power  necessary  to  perform  these  tasks  with  any  degree  of  an¬ 
alytic  computation  is  far  higher  than  that  available  on  a  typical  desktop  computer,  and 
so  a  purely  numerical  solution  would  be  much  easier  to  compute,  if  less  accurate  than  an 
analytic  solution.  This  involves  being  able  to  take  Laplace  transforms  numerically  as  well 
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as  improving  upon  the  accuracy  of  the  numeric  computations  already  implemented  in  the 
code  provided  in  Appendix  A. 

5.2  Conclusions 

Based  on  the  analysis  of  the  modified  wave  equation  performed  in  Chapter  III,  it 
is  possible  to  filter  acoustic  signals  analytically,  yielding  a  transfer  function  that  can  be 
applied  to  the  Laplace  transform  of  any  signal.  With  this  altered  signal,  the  original 
signal’s  Laplace  transform  can  be  found,  which  can,  in  turn  be  inverted  to  yield  the  exact 
original  signal. 

By  generalizing  the  number  of  signals  and  receivers  in  a  room,  the  problem  of  filtering 
acoustic  signals  can  also  be  solved  exactly,  given  that  certain  conditions  are  met.  These 
conditions  relate  directly  to  the  determinant  of  the  matrix  of  transfer  functions.  If  the 
sources  are  too  close  to  one  another  or  the  receivers  are  too  close  to  one  another,  then 
the  determinant  of  the  matrix  is  close  to  zero,  making  the  matrix  ill-conditioned.  If  the 
matrix  is  not  ill-conditioned,  however,  the  computation  of  its  inverse  is  stable  and  the 
unique  solution  of  a  vector  of  the  original  signals  can  be  produced. 

Further  generalizing  the  problem  of  finding  and  inverting  the  room  transfer  functions 
to  accommodate  mobile  sources  results  in  a  finite  convolution  with  more  terms  than  with 
a  non-mobile  source,  which  in  turn  requires  taking  two  Laplace  transforms  to  result  in 
an  equation  for  the  original  signal.  This  doubly  transformed  signal  can  be  transformed 
by  taking  the  inverse  Laplace  transform  with  respect  to  the  second  transform  variable, 
and  then  with  respect  to  the  first  transform  variable.  The  analytic  form  of  this  signal 
is  very  complex,  if  it  exists  at  all.  The  inverse  Laplace  transform  may  have  to  be  taken 
numerically  in  order  to  yield  a  result. 

The  difficulty  in  implementing  these  analytic  solutions  is  quite  clear  from  the  graphs 
of  the  transformed  received  signals  presented  in  Chapter  IV.  Because  of  the  enormous 
number  of  calculations  required  to  compute  the  transfer  functions  of  the  room,  a  reliable 
solution  cannot  be  found  using  the  code  in  Appendix  A,  and  so  a  better  program  must  be 
developed  for  finding  these  transfer  functions  numerically. 
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5.3  Recommendations  for  Future  Research 

Given  the  sheer  volume  of  assumptions  made  in  this  thesis,  there  are  countless  av¬ 
enues  of  research  opportunity  available  in  this  field.  The  assumption  of  point  sources 
and  receivers  is  highly  invalid  in  real-life  acoustics  when  dealing  with  human  vocal  signals. 
Since  the  human  head  is  not  best  described  as  a  point  source,  the  signals  propagated  from 
a  human  mouth  are  not  omnidirectional.  The  shape  of  the  throat  and  mouth  acts  as  an 
antenna,  shaping  and  directing  the  signal  in  the  direction  the  head  is  facing.  Furthermore, 
microphones  are  of  a  definite  size,  and  so  they  cannot  be  considered  as  point  receivers,  but 
as  directionalized  receivers.  After  all,  human  ears  act  as  antennae  in  receiving  the  signals, 
making  the  reception  directionalized.  Modifying  the  wave  equation  in  a  way  that  accounts 
for  these  changes  in  signal  sources  and  receivers  would  yield  a  more  accurate  result. 

Exploring  the  nature  of  the  matrix  of  transfer  functions  is  another  area  of  possible 
research  opportunity.  If  a  matrix  is  ill-conditioned,  or  is  even  noninvertible,  other  methods 
than  inversion  are  possible  to  find  the  solutions.  Using  the  Moore-Penrose  pseudoinverse 
is  discussed  in  Chapter  III  as  a  way  of  finding  these  solutions,  but  if  the  matrix  is  ill- 
conditioned,  even  these  pseudoinverses  require  other  transformations,  such  as  a  singular- 
value  decomposition  or  a  QR  decomposition.  Finding  a  method  of  solving  problems  of  this 
type  could  be  crucial  in  situations  where  the  location  of  the  sources  and  receivers  cannot 
be  changed,  but  must  be  accepted,  such  as  the  intelligence  gathering  scenario. 

More  importantly,  the  research  in  the  aforementioned  fields  is  useless  without  a 
proper  implementation.  A  more  efficient  computer  program  must  be  developed  based  on 
the  findings  of  whatever  research  is  conducted.  Perhaps  using  a  different  programming 
language  would  be  beneficial  to  the  efficiency  of  the  method,  such  as  using  C,  where 
the  programs  can  be  compiled  to  run  faster.  Once  the  programs  are  running  faster,  it 
is  easier  to  implement  more  accurate  computations,  where  the  transfer  function  can  be 
computed  with  more  than  two  million  cycles,  as  it  is  done  in  the  code  written  for  this 
thesis.  Furthermore,  the  use  of  experimental  data  will  be  very  helpful  in  implementing 
these  codes.  Not  only  would  the  implementation  be  used  to  actually  invert  the  effects  of 
the  room  on  a  signal,  but  the  actual  data  would  include  the  directionalized  nature  of  both 
the  sources  and  the  receivers,  giving  a  more  accurate  result.  The  lack  of  experimental 
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data  in  this  thesis  greatly  diminishes  the  applicable  value,  but  the  groundwork  for  future 
experiments  is  now  laid. 
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Appendix  A.  MATLAB  Code 

This  appendix  contains  the  code  for  finding  the  transfer  function  of  a  room  and  produces 
the  signal  received  by  a  microphone  located  at  the  user  specified  location..  The  first 
section  contains  the  program  for  a  single  source  and  single  receiver.  The  second  code  is  a 
program  that  finds  the  Laplace  transform  of  a  signal  numerically. 


A.l  Single  Source,  Single  Receiver  Transfer  Function  Program 


clear 

clear 


H=3 ; 


L=10 ; 
W=10 ; 


xrl=2 

yrl=2 

zrl=2 

xsl=0 

ysl=9 

zsl=2 

bx=10 

by=10 

bz=10 

c=374 


9; 


v=-20000: 10: 20000; 

for  i=l : 4001 
sum(i)=0; 
sumplus (i)=0 ; 

end 
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for  i=l : 4001 


for  1=1:10 

for  m=l : 10 

for  n=l : 10 

num= (sin(n*pi*xrl/W) +(n*pi/bx/W) *cos (n*pi*xrl/W) ) * 
(sin(m*pi*yrl/L) + (m*pi/by/L) *cos (m*pi*yrl/L) ) * 
(sin(l*pi*zrl/H)+(l*pi/H/bz) *cos (l*pi*zrl/H) ) * 
(sin(n*pi*xsl/W) + (n*pi/W/bx) *cos (n*pi*xsl/W) ) * 
(sin(m*pi*ysl/L) + (m*pi/L/by) *cos (m*pi*ysl/L) ) * 
(sin(l*pi*zsl/H)+(l*pi/bz/H) *cos(l*pi*zsl/H)) ; 
den(i)=( ( ( . 5+bx"2*W~3/ (2* (n*pi) "2) ) * 

( . 5+by"2*L~3/ (2* (m*pi) ~2) ) * 

( . 5+bz~2*H~3/ (2* (l*pi) ~2) ) ) * 

(v(i)"2+c~2*(  (l*pi/H)  ~2+(m*pi/L)  ~2+(n*pi/W)  ~2) ) )  ; 
sumplus(i)=sum(i)+num/den(i) ; 
sum(i)=sumplus (i) ; 

end 

end 

end 

end 

A.  2  Numerical  Laplace  Transform  Program 

clear 

t=0 : .0001:5; 

for  i=l: 50001 

chirp (i)=sin(2*pi*t(i) ~2 ) ; 

end 

for  i=l : 4001 

lchirp(i)=0; 
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Ichirpplus(i)=0; 


end 

v=-20000: 10: 20000; 

for  i=l : 4001 

for  j=l : 50001 

lchirpplus (i)=lchirp (i) +exp (-v(i) *t  ( j ) ) *chirp ( j ) ; 
lchirp(i)=lchirpplus (i) ; 

end 

end 

figure 

plot (t , chirp) 
figure 

plot (v, lchirp) 
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